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Temperature Distributions within the Earth.* 


By 


Haruo Mixt 


Geophysical Institute, Faculty of Science, Kyoto University, Kyoto, Japan. 


Abstract 


Several features of the temperature distribution within the earth’s interior were 


derived from the modern theory of solids. 


The main results are as follows: (1) 


GRUNEISEN’s parameter in the B-layer (33-413 km.) must be greater than 1.5. (2) Tem- 
perature in the D-layer (1000-2898 km.) increases almost linearly with increasing depth 
and its gradient depends on the mean atomic weight rather than on GRUNBISEN’s 


parameter. 


For instance, the gradient is 1.79°/km. for A=20 and 3.74°/km. for A=40. 


(4) The temperature gradient in the B-layer decreases with increasing depth. 


Notations. P: Pressure, 7: Distance from 
the earth’s center, g: Acceleration of gravity, 
po: Density, T: Temperature, Kr: Isothermal 
bulk modulus, Ks: Adiabatic bulk modulus, 
a: Coefficient of thermal expansion, Cp: 
Specific heat at constant pressure, C,;: Specific 
heat at constant volume, ¢=Ks/p, A: Mean 
atomic weight, V: Molar volume, 7¢: Griin- 
EISEN’S parameter, v,: DEBYE’s maximum 
frequency, N;: Loscumipt’s number, vp: Velo- 
city of longitudinal elastic wave, vs: Velocity 
of transversal elastic wave. 

Theory If we assume that the condition 
‘of hydrostatic equilibrium is satisfied within 
the earth, that is, dP=—gpdr, then the varia- 
tion in density with depth is expressed by the 
following equation 


do_ yet egos 
dr \@P/7 adr OT /p dr 
ae Gale 
© fae Phas (1) 


Let the radial temperature gradient be 
expressed as follows, 
aT _Ta dP _ 

dr ~ pCp dr 
where ct denotes the difference between the 
actual gradient of temperature and the adia- 
batic gradient, —Tag/Cr. In virtue of the 


Fee = (2) 


thermodynamical relation between isothermal 
and adiabatic bulk moduli, 


Kr7/Ks=1—Ta?K7/0Cp, 
the equation for the density variation with 
depth may be written as follows (Brrcu, 1952): 
do _ go” ae 
dr Kr oCp 


heer Fe, FoF 
(3) 


Let the earth be divided into thin spherical 
shells, the inner and outer radii being 7 and 
7. We use the suffixes 1 and 2 also to dis- 
tinguish the physical quantities of the ma- 
terials existing at 7, and 7,. Then (3) can 
be approximated by the following equation 
(Mixt, 1952): 

P1— 0s J Pit Pe Oi +02 

it ee + at— fo 

Y iad i} 7) 2, 2 
which can be transformed as follows: 


=—(7 5) 
UWE} Va Vz 
1 gfA A 1 € a) 
= — - (a4 ; 
2 ee gs AAG ie 
Or 
1 4S 1 ( alt a) 
L— ——| | an : 
Tlie Ts ( Vz 2 is i o) 


Therefore, we get 


’ 
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The ratio of molar volumes V,/V. can be 
derived from the definition of GriiNEISEN’S 
parameter, thus 


(9) 


the temperature can be expressed as follows: 


yi 2! ih ) 
je a (8) 
= 0 ter 


We have known the value of ar at each 
depth, that is, we have known a’s as a func- 
tion of t at various depths. From this func- 


- 2 endl 
a 9 


Dis ~ Use tional relation and equation (8), the corre- 4 
sponding values of T and ct can be known. 
ats In other words, we can know 7t as a function 
d log vin 9 N, 1 + of Tf. Thus we can find the distribution of 
rae Vv? and m=(Foy 1 9 ) temperature in the earth by numerical inte- 
ee bare gration of equation (2) after substituting the 


Since vp and vs are quantities which can 
be known from observations of seismic waves, 
V,/V, can be calculated from (5) taking re 
as a parameter. Then ar can be estimated 
from equation (4), since @ is a known quantity 
and g is nearly constant within the earth’s 
mantle (about 10° cm. sec.~?). 

On the other hand, in virtue of the follow- 
ing relations: 


Results 


above functional relation between T and vt. 
The first term of the right-hand-side in equa- 
tion (2) is the adiabatic temperature gradient 
and it can be transformed into — eae by 
(7), and the coefficient of ZT can be calculated 
from experimental values. C, at high tem- 
perature is 2.49432 x10® erg. deg.-! mol.-! by 
Dutonc-Petit’s law. The results of calcula- 
tions are given in Table I. 


Table I Absolute Temperatures within the Earth 


Mean atomic weight 20 


Mean atomic weight 40 


(km.) Veemelico) V¥ac2.0 Ve — 1.0 ¥¢ =2.0 
33 1000* 1000* 1000* 1000* 
B 100 982 1058 960 1122 
200 9338 11389 855 1292 
800 819 1190 616 1397 
— 413 
a) 
—— 1000 2000%* 2500* 8000%* 2000* 2500* 3000* 2000* 2500* 3000* 2000* 2500* 3000* 
1200 2256 2752 8248 2340 2825. 3318 2525 3021 3516 2721 3210 3698 
1400 2590 3081 3572 2714 3188 3669 3211 3701 4191 8511 3988 4464 
1600 2952 3488 3923 3097 3560 4030 8954 4489 4923 4318 4782 5247 
1800 3348 3828 4308 3492 3944 4404 4766 5244 5722 5147 5601 6055 
D 2000 3711 4187 4662 | 8868 4311 4761 5507 5981 6454 5933 6376 6820 
2200 4136 4606 5076 4270 4708 6148 6376 6843 7311 6770 7204 7637 
2400 4439 4905 5372 4616 5041 5473 6990 7454 7918 7487 7912 8338 
2600 4813 5275 5737 4988 5405 5829 7751 8210 8669 8257 8674 9091 
ae 5198 5656 6113 5322 5770 6186 8533 8987 9442 9028 9487 9846 
: | | 
a a 


K 


assumed values 


Temperature Distributions within the Earth. 


C-layer (413-1000 km.) was not considered 
in the above calculations. 

If GriinEISEN’s parameter is taken as 1.5, 
the temperature decreases with increasing 
depth in the B-layer. This in an inconceivable 
conclusion. Therefore the value of GriiNEISEN’s 
parameter must be greater than 1.5 in the 
B-layer. 

We see that the temperature increases ap- 
proximately linearly with increasing depth in 
the D-layer. The mean temperature gradient 
in this layer is as in Table II. 


Table Il. Temperature gradient in 
the D-layer (°/km.) 


Mean Ye =1.5 
atomic Total 


weight. | “esr 000 = 2000 2500|3000 2000 2500 30 


| 1 =2.0 
| | 


00 


20 
40 


Table II shows that the temperature gradi- 


ent depends mainly on the mean atomic 
weight rather than on GriinEISEN’s parameter. 

The temperature gradient in the B-layer is 
as in Table III. 


Table HI. Temperature gradient 
in the B-layer (°/km.) 
2 Mean atomic Mean atomic 
Dept (kins) weight 20 weight 40 
0.864 1.761 
200 0.661 1.375 
300 OL Sail 0.638 


Table III shows that the temperature gradi- 
ent in the B-layer decreases with increasing 
depth. 
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On a Self-Exciting Process in Magneto-Hydrodynamics (III). 


By 


Hitoshi Takrucar and Yasuo Saimazu. 


Geophysical Institute, Faculty of Science, Tokyo University, Tokyo. 


Abstract 


In our recent papers (referred to as Paper I and II), it has been shown that a 


self-exciting process is possible by which 
produced and maintained. 


the earth’s main magnetic field may be 


In these papers, however, a kind of theoretical “cut off” has been made, and 
only the magnetic fields of comparatively lower harmonics have been taken into acount 
in the discussions. In the present paper, the relative importance of the magnetic 


fields of higher harmonics in driving the self-exciting dynamo is estimated. 


It is shown 


that the results in Paper I and IT remain true even if we include the magnetic fields 


of higher harmonics in our dynamo model. 


§1. In a recent paper, E.C. Butuarp and 
H. GeLLMAN (1952) have made comments on 
our first paper (1952) on a self-exciting process 
in magneto-hydrodynamics (which will be 
referred to hereafter as Paper I) and stated 
as follows: “The existence of such eigen- 
values does not establish the existence of 
dynamo solutions of Maxwell’s equations, 
since in Figure 1 (in our Paper I) all har- 
monics beyond the second have been omitted 
....though naturally the existence of solutions 
of the simpler equations makes it more likely 
that the complete set also has solutions”. In 
‘order to make the matter clearer, we shall 
reconsider the way of reasoning taken in 
Paper I. In Paper I, the existence of velocity 
fields (relative to the earth’s mantle) of T,°- 
and S.2°-types in the earth’s core was assumed 
a priori. Owing to the inductive action due 
to these velocity field, groups of magnetic 
fields may be caused and maintained in the 
earth’s core. In these groups, there may 
exist one which contains a dipole (S,°-type) 
magnetic field as its constituent. It is this 
group that we studied in details in Paper I 
and in our more recent paper (1952) (which 
will be referred to as Paper II). In the model 
earth adopted in our study, the earth’s core 
is assumed to be surrounded by an insulating 


mantle. In the mantle, only the magnetic 


field of dipole type is assumed to exist. In 
Paper I, we investigated under what condition 
the group of magnetic fields considered above 
can exist. The problem is reduced to an 
eigen-value problem for « (conductivity of 
the earth’s core) times V (absolute value of 
the velocity). The eigen-value problem is 
solved in Paper I and II. In view of the 
results obtained in these papers, we may ex- 
pect that the group of magnetic fields con- 
sidered above can exist in the earth’s core. 
In the course of the analyses, however, we 
have made a theoretical “cut off”. BuLLarp 
and GELLMAN’s remark above referred to is 
on this point, and may be stated as follows. 

There are infinite numbers of magnetic 
fields which are coupled with the dipole mag- 
netic field by the inductive action of the 
above fluid motions. These fields are labelled 
with the suffixes and m, which are the 
degree and order of the surface spherical 
harmonics belonging to the fields respectively. 
These fields are caused and maintained in a 
self-exciting way. From the physical point 
of view, the self-exciting dynamo may be 
considered to be driven by the magnetic fields 
of comparatively lower harmonics. In fact, 
in our Paper I, only the magnetic fields of 
n, m<2 were taken into account. It is, how- 
ever, not self-evident from the beginning 


6 


that the magnetic fields of higher harmonics 
play no important roles in the self-exciting 
dynamo. In this meaning, BULLARD and 
GrLLMAN’s remark is to the point. It is the 
object of the present paper to estimate the 
relative importance of the magnetic fields of 


S-type (Poloidal type) magnetic field 
—n(nN+1)Sr™(r)r271¥n™ , 
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higher harmonics in driving our self-exciting 
dynamo. 

§2. In order to investigate the general 
theory of magneto-hydrodynamical couplings, 
we shall make use of the following expressions 
for the (7, 0, ¢)-components of magnetic fields. 


DSn™ m Ign —- ,OY¥n™ 
fe -|r Sn" + (n+V)Se i¢ fa 
| 1 OYn m 
— w ie 
sin 60d” 
0, 
@S,™ dSn Oe 
Paces » fae rn ae £m 
curl Hs= lr 7 as ae “} sin 00¢ ’ 
a ; ae 
00” 
—nin-+) Oo ue EY) ae Yn™, 
Y dr 
oTy. as GeSy a Sa OY Gn 
Vd g=+-= TE 3 1 2n j= LR A) Vs (i Waa 
eee eae drt hope? 
OG 
ai r m—1 = 
| fr he (2.1) 
T-type (Toroidal type) magnetic field 
0, 
' OY. m 
—T.™ yr” n 
H7= ee sin 00¢ ’ 
\ T,™n)r° en : 
ee Va m : 
ra i 34) tes OV Re 
curl Hr= ie ue fr an 
; n—-1 0Yn™ 
uw 13 5 
sin 00g ’ 
a 
Pane: Lee ene) ed LOLne 
r"— 
dr’ r dr sin 60¢ ’ 
wolf Y ve Oya 
ad.” 8 ee ee ae (222) 


The velocity field in by earth’s core is assumed, as in Papers I and II, to be 
an arbitrary linear combination of the following two fields. 


S-type velocity field 
—Mn+1l)V okt), 
adv; Me m 
Vows -[r oP +(n+)Vsh fp a 


=| io 


sin 0p’ 


expressed as 


Yn” 


LOX n 


see eee 
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with m=2 and m=2c (2c means cos 2¢ in the 
expression of (2.5). 


T-type velocity field 


0, 
ia m "(r) n OY ,™ 
Vz= sin 00¢ ’ 
joe n™ 
Va,n(“)r ie) (2.4) 


with m=1 and m=0. In (2.1)-(2.4), Yn™ is a 
surface spherical harmonic with degree » and 
order mm, 


vudvU 


Yn" =Pym(cos 8) ind (2.5) 
The expressions for the magneto-hydro- 
dynamical couplings of general (poloidal and 
toroidal) magnetic fields with general (poloidal 
and toroidal) velocity fields will be given as 
follows. 
With the suffixes 
@: primary (inducing) magnetic field, 
B: velocity field, 
ry: secondary (induced) magnetic field, 
the orthogonal condition for the Maxwentu’s 
equation ((3.6) in Paper I) is expressed by 


| | | SIp?Ha+4ie curl (Ve x Hx)]-Hyr sin 0drdddg=0 . (2.6) 


Executing the integrations with respect to @ and ¢, we get 


CL clL——“aQ 


where 
na 2n 
— F 
Saoy — | | 
my) 
0 vu 
te Qn 
ee 

IDS 

0 0 


[Su55 25° S,17407=0 : 
a a 
[SSeTy4+ ST .Tydr=0, (2.7) 


[7° Hsa + 4zx« curl (Si) Vex Hsa)]-Hsy sin 0d6d¢ , 
6 


[7?Hra + 4r« curl (3) Vex Hra)|-Hsy sin 0ddd¢ , etc. (2.8) 
B 


a, LoS ea ame P . . : 
We shall call SaSy, TaSy, SeTy and TT, the coupling coefficients. Changing the independent 
and dependent variables by the following equations 


7 —Ces (a: radius of the earth’s core), 


(ean ~ il (e® ie 


Vg"(7) a aN Vig n, mV gin n (E) 


Tair Tum) 
Teel E ). (2.9) 
VO)! O\Vr,nmVrME) 


~ we get the following expressions for the coupling coefficients. 


Couplings by S-type velocity field: 


SeSy=[Ma@, 7)na(Mat fs yK (ar, B)ne(np+ gi —yKUrB, nna a+ Vg _|ny(My+ 1)SyE"~! 
AG, Y2 ) Ny —2 
+| He 1) Na( Na+ 1) eae ) +yK(ar,B)ne(ne+1)- ae oe —yK(7B,Q)Ma(Ma +1)- (926° | Sen, 


dé 
(2.10) 
1 LD)SyEm=14 yL (rar, B) np (np-+ 1) LE? Sem, (2.11) 
TaSy =yL 7a, B)ne(mp+1)gsny(My+1)SyE"Y EY B de 
ae tL ioe -9 YI 2) 
S.7y=| bia, +) Ny(My+1)SaVsp&"et"p-? +yL(ar, B)ne(mp+1) 


+yL(7B, O)na(Ma+1) a(928") foe, (2.12) 


dé 
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FaTy=| Het rymalne + Df YROB, r)rylor-t Ta Vana e 
+yKar, Byung +) oe > ge. (2.13) 
Couplings by T-type velocity field: 
S,Sy=[Ia, p)na(nat fs tel (Br, O)na(na+ galery +1)SyE"y? 
+| Na. rytie(na+ ILE) 4 eX(B7, Matta + 1) 2 ae ise 2 (2.14) 


dé d& 
T= (2.15) 
7 : = na 1) HOE”) pygmy (2.16) 
S.Fy=—| 2K, r km Sa roe "9 +2K(7B,Q)na(nat 1) dE ie tess { 
Tu Ty=[Hax, 1)na(ta+l)fr+2L(aB, 7)iy(My+1)T a VrpEra tp] TE". (2.17) 
where 
y=4rKaVy, oc, v=ArKraV,,o, S= ae +(n+1)S, 
= ed Vg Pe 241) A sas 
a Cn Es ae 
91 = Sa VspE"a*"p-3, 92=SaVspE"a*p-3, 93= Ta VspE"«*"p-2 
Gea a Vierpbtorne et, i (2.18) 
and 
I(a, n=\\ YaYy,sin 6d6éd¢ , (2.19) 
00 
nw 2 
(2422 if. OY OY se )v eee an 
aed ee Sn SUG SOD oe Ae 
00 
ile (e3 0¥e OY. ee 
Y Ad0dé . 2.21 
Sat hy ahi @0 sinéd¢d sinddg? oes) * meno? 22) 
uv 0 


In (2.19)-(2.21), Ya, for example, denotes Y,". ae Taihaihadine To reat S.-types Gi. 


cities and in which S,°-type magnetic fields is 
found as an element. In Fig. 1, is shown 
the chain of couplings of magnetic field up 


The selection rules for the coupling coef- 
ficient integrals in (2.19)-(2.21) are as follows. 
As to the integral J(@,7), we must have 


Na=Ny, Ma=My. (2.22) to the fourth harmonics terms. It is needless 
Otherwise the integral vanishes. A(@f§,7) to say that the chain is continued infinitely 
vanishes unless with infinite numbers of magnetic fields. 


Natne+my=even, mackmp=my. (2.23) As is seen in Fig. 1, magnetic fields of S,°, 
T.°, T° and T.°8 types make a closed circuit 
of couplings with harmonics up to the second 
degree. This simple chain was studied in 
These are the selection rules for the magneto- Papers I and II. We shall call this chain 
hydrodynamical couplings. We shallnow study the first approximation chain of the self- 
the magneto-hydrodynamical couplings, which exciting process. Magnetic fields of S,°—, S,2°- 


L(aB, 7) vanishes unless 


Nat+np+ny=odd, matkmp=my. (2.24) 


=) 
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Fig. 1. Chain of induction mechanism by 7, and S,2°-types fluid motions. 


and S,25-types become the elements of the 
coupling chain if the harmonics up to the 
third are taken into consideration. This step 
may be called the second approximation. 
Including harmonics up to the fourth, the 
chain expands and includes twelve types of 
magnetic fields. It will be called the third 
approximation chain. The numbers of ele- 
ments of coupling chain increase rapidly as 
the higher harmonics terms are included. 

§ 3. Since the coupling chain is much com- 
plicated, it is hopeless to solve the Maxwell’s 
equation with higher harmonics terms in 


in which the T,°-type velocity field tends to 
infinity in magnitude. As we saw in Paper 
II, this special case is of some geophysical 
interests. This corresponds to the case x00 
in Paper I, and we got the eigen-value Y for 
this case in Papers I and II, taking only the 
harmonics up to the second into consideration. 
We shall now try to get the eigen-value Y 
for this case, taking the harmonics up to the 
fourth into consideration. The condition zoo 
simplifies the treatment of the higher ap- 
proximation considerably. This circumstance 
will be shown below for the first approxima- 
tion already treated. From (5.2) in Paper I, 


general. We shall consider the special case 
we have 
et IVT 
E a See =< ae 5940 ope VOT 2) + 5 & VET | 
gO TE ATE SA) CE VOT + SEVIS |-2eVKET , 
a ans = 2a Vi9ET +2 p36 Vi val a + Porte {e @veset | 


(3. 


At first, we shall assume that S,° magnetic field is of finite order. From the first equation 
in (3.1), we see that T,* field is of finite order. Making aco in the fourth equation, we 
see that T,2° field becomes zero. T,° field becomes the order of x by the second and also 
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by the third. In short, we get T,°~x and 
T.2°=0 when «>. Taking the similar 
way in the third approximation, we get Tot, 
Ti ~e and T°, S,°+T 29, Tf°-and T° =0 
Thus, the twelve simultaneous 
equations are reduced to seven in the third 


as woo, 


shown in (3.3) below. 
Putting 
fears rn 
T S=2T (3.2) 
we have T,° and T,° of the finite order. 
Making use of (3.2) and the results in SZ. 


BEDE rune uon: The seven cme will be we have 
ee 4 ie COS em 
ae ERs 5 
486%? yf (480 54520) vy04.s(ES TE + 3¥°) [P=0, 
7 dé dé 
Ge ,@ Var ol ee 6 at") 48 x6 ¢ stVi ee 0, 
Sige de Ede ae aE 
288K6 | = o/-dVi® )+ a {7 a0Re =| 96 X6 7H os 
ee air 3V2°) +2 (FV, Te 
oo a Geers OVE he ; 
96 x 12 ee. 32 x6 E ‘i a 3V,°) 2,4 (T 9 Vo2°E*) |= 
288 x 6 y 5 , 48 ae 8 dS;° ) 
OE 2S VY, 2 ee = = 
my TRV See te aE 
eo, ee “-+452°) Vi —25e(E or +30") |- 480 x 6. yT 28 VE =0 
7 dé dé 
EBX yF V2 — 0x1 5 20,8 — ee y yT OV 2°H=0, 
p20. a" 1. 80/ OTe dT, 
Sea (ES SZ a) 
Pies Gee dE iit ce 
_ 192x5 aT es oh Gay 36 _ 480 x 20 gcd 
: uf ar o( EST +3V,1)—38 vet V, ee eS 
160 x 18 wire ava oH \ ea d 20 a 
a — V.2° \E2+-17—(T 1° V.2°E*)- 160 x 207/78 VPS=0 
Zs Tutt Ven 430s) 417 S(T eVect), [4 LVN 


These are the simultaneous differential equa- 
tions for the radial functions. 

Assuming that all magnetic fields in our 
earth’s mantle vanished at roo, we get the 
following boundary conditions 


v2.0 dV2E) 
V9 j= V. ( =() 

(E)= Vi2(E) = dé 
dS mm 

Sp” ‘sg E i 

(&)4 Ont1 dé ; 
If ni™E) = 
at g=1, (3.4) 


Thus our problem is reduced to solve (3re)) 
under the boundary conditions (3.4). It is, 


see eee 


however, almost hopeless to solve this eigen- 
value problem for y=4z«aV2,.0 by the trial 
and error method. We shall solve it in an 
approximate way, which is practically the 
same as that in Paper I. 

At first, we shall replace the equations in 
(3.3) by their equivalents of the types in (2.7). 


As in Paper I, we shall put rather tentatively 
ViX@)=1-E 5 VS) =A =e) 
Next, taking the boundary conditions in 


(3.4) into consideration, we take the following 
trial coordinate functions. 


Si(6)=4,°(1— ao , 
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BACEy WNP Ai 
emis ey 
F.%(E) By 
T.8(E)\_{ By \ 
T &é) : BY (8). 
Te)! \ Bes 


(3.6) 


where A,°, .... are undetermined constants. 
Since only one adjustable parameter A,” (or 
B,™) is given for each coordinate function, 
the present approximation using coordinate 
functions in (3.6) may be called the first ap- 
proximation in the sense used in Paper I (see 
page 7 in Paper I). 

The approximation using A,°, A;°, A;2°, 
B,°, B28 is called, however, the second ap- 
proximation in the sense used in §2 of the 
present paper. Similarly, the approximation 
using A;", A,, A,**, B.°, B25, BY, BZ is 
called the third approximation. In any case, 


Table I. Equations for 


Magneto-Hydrodynamiecs (III). 11 
inserting (3.6) into the equations of the types 
as in (2.7), we get the equations of the fol- 
lowing types 


Sy eee 96 
35 35.x7 


All of these equations for the third approxi- 
mation (in the sense of the present paper) 
is shown in Table I. This table is a generali- 
zation of Table I in Paper I. In order that 
these equations are compatible with one an- 
other, the determinant formed by the coef- 
ficients must be equal to zero. By putting 
the determinants at the first, second and the 
third stages (in the sense used in the present 
paper) equal to zero, we can determine the 
eigen-value 


YB." + yA7O=0, etc. 8:2) 


y=4rKaVy,20 (3.8) 


for each stage. The values of Y thus calcu- 
lated are shown in Table II. 


the third approximation. 


A,o Bo ByS A,o A326 Bo | Bps 
F ‘ : 24 | 96 0 
=? : 35" : 35.77 : 
o = 0 ce i 0 0 0 
30 35 35.55 
72 24 _ 48 48 ' 
v 35.5B 35.5 0 35 11 36.11.13” 
“eo oenige | 24.8 | 320 0 _ 24.12, 
v e Sa Bp” ~~ 49 77.13 4 49.13” 
ci ‘ Semel | 24.8 96 i 
? 35.11 7 v 8 ru 49.13” 
; Eee. 80 
| 0 0 
: B q 11.13 | 99 
=. i 800 


Eigen-value Y by the 
approximate method. 


Table IT. 


Ww 
Ist 68.51 
2nd V7 .79 
8rd 80.91 


The value Y=68.51 is nothing but the one 
obtained in Paper I. Other values of Y in 
Table Il do not differ much from this. We 
may infer from this result that the self- 
exciting dynamo is driven by the fields of 
comparatively lower harmonics. Thus atte 
results in Papers I and II are shown to remain 
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true even if we include magnetic fields of hydrodynamics”. Journ. Phys. Earth, 


higher harmonics in our dynamo model. This ib, (Gis); 
result may be taken as an answer to Butnarp TAKEUCHI, H. & SHimazu, Y. 


and GELLMAN’s comment referred to at the egoe oe eke ee eee ee 
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1952 “On a self-exciting process in magneto- 
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Convective Fluid Motions in a Rotating Sphere. 


By 


Hitoshi Takeucut and Yasuo SHimazu. 


Geophysical Institute, Faculty of Science, Tokyo University, Tokyo. 


Abstract 


In connection with the self-exciting dynamo problem studies in our previous papers 
(referred to as Paper I, II and III), the present paper deals with the stationary con- 


vection motion of fluid within a rotating sphere. 
In §5, 6 and 7, the problem is reduced to 


are made for plane boundary problems. 


In §3 and 4, preliminary surveys 


solve the equations (6.16)-(6.18) under the conditions (7.15) and (7.16). This is an eigen- 


value problem for the Rayleigh number {= 


eigen-value problem is solved in §8 in an approximate way. 
the fluid motion required in Paper I-III is possible in the case when 2~102. 


atabg 
key 


ah wa ; 
containing a parameter 2=——. This 
Vv 


In §9, it is shown that 
The 


corresponding kinematical viscosity » of the earth’s core comes out to be about »y~10 


poise. 


§1. In previous papers (H. Takrucnr and 
Y. Saimazu, 1952, 1953), which will be referred 
to hereafter as Paper I, II and III, the possi- 
bility of a self-exciting process was established 
by which the earth’s main magnetic field is 
produced and maintained. The self-exciting 
process is considered to be maintained by the 
induction currents caused by the motions of 
the fluid of which the earth’s core is com- 
In order to make the study on the 
self-exciting dynamo complete, however, there 
remains another problem to be solved. That 
is the problem on the fluid motion itself. Can 
the required fluid motion really take place in 
Under what condition is 
that fluid motion possible? Is the admissible 
fluid motion appropriate to maintain the self- 
exciting dynamo under question? These 
questions will be considered in the present 
paper. 

§2. All possibilities having been considered, 
it is now believed that the fluid motion within 
the earth’s core is that of convection due to 
the non-homogeneous heating to which it is 
subjected. (E. C. Buniarp, 1949; W. M. 
Exsasser, 1950). Our previous studies in 
Paper I-III are made on this convection cur- 


rent model. The S.2°-type velocity field 


studied in these papers is nothing but the 
mathematical expression for the convection 
current. It is this S,?°-type velocity field that 
makes our self-exciting dynamo possible. In 
view of these circumstances, our immediate 
problem to solve may be stated as follows. 

(1) Can the fluid motion of S,?°-type exist 
stationarily within the earth’s core? If it 
is, under what condition ? 

(2) Also it must be shown that under the 
condition studied in (1), the S,?°-type velocity 
field is the easiest mode of motion to be ex- 
cited. Otherwise, the S.?°-type velocity field 
may be overcome by fluid motions of other 
more powerful modes. 

The stationary velocity field of S,?°-type 
reminds us of the well-known Benarp-cell. 
It may be that the fluid motion in the earth’s 
core is essentially of the similar character as 
the Brnarp-cell motion in a vessel which is 
heated from below. This will be the main 
point of view held in the present study. Many 
studies have been made about the convective 
motion which a viscous fluid undergoes when it 
is heated from below. In most of the studies, 
the problem is treated for plane boundaries, 
and the effects of earth’s rotation upon the 
fluid motion are not taken into account. The 
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effects of magnetic fields are also neglected. 
In the next two sections, we shall take the 
earth’s rotation and the existence of magnetic 
field into account, and treat the problem for 
plane boundaries. 

$3. We shall consider the fluid of density 
0, viscosity y, kinematic viscosity ame 
thermal diffusivity k and coefficient of thermal 
expansion a. The fluid will be considered as 
being contained in two parallel walls placed 
horizontally at a vertical distance of h. The 
fluid is heated from below. The temperature 
gradient in a state of no convection is BP. 
Referring to the rotating rectangular coordi- 
nates as is shown in Fig. 1, we have the 
following equations. 


Fig. 1 
Du_ ‘Op ; 
0 ie oy + pp-ut2owv , 
Derr. OD 
ODE ay + up?v—2owu , 
Dw Op 5 
(Disumelie a) ee 
Dae : 
Dt Sarl . 


Customary notations are used. Attaching the 

suffixes 0 to the quantities in the states of no 

convection, we have 
M=VUp=Wo=9, Po= Plz), 
aT) _ 


To = T)(2), 


ae 2, Po= 00,001—aT)), 
dpy 
——~—pg=0 
vias ale (3.2) 


?o,9 being a certain mean density correspond- 
ing to a certain mean temperatures (or 2). 
We shall put 


p=poth, T=T.+T’, p=ate’, (3.3) 
for p) £ and p im (3.1); and neglect the se-, 
cond order terms in w, v, w, p’, T’ and 0’. 
Furthermore, we shall assume that the density 
variation of the fluid is caused by the thermal 
expansion alone and the fluid behaves other- 
wise as if it were an incompressible one. Thus 


we get 
0=0,0(1—aT), op =—000.AT, (3.4) 
Ou Ov Ow 
Ou, Ov Uy 3.5 
Ge By * Oz 2 


Inserting (3.2)-(3.4) into (3.1) and denoting p’, 
T’ and ~o,9 as ~, T and p anew, we have for 
the state of stationary convection 

at 2? +-20v=0 ; 


ve 


ypu 


ype Ope =) : 
0 Oy 


See £ 
ypy7w A qe ee 
kp?T+pw=0. (3.6) 


While the above equations are derived for an 
incompressible fluid, H. Jerrreys (1930) has 
shown that the same equations can also be 
applied to compressible fluid provided the 
density does not vary greatly within the 
system and provided also we interpret 8 as 
the temperature gradient in excess of the 
adiabatic one. Putting 


CE OE, 7a, 109) 


£ (many) 
= e” 


(ule), v2), we), D@),. 12), ae 


we get five ordinary differential equations for 
five unknown functions w(z), v(z), ...Putting 


mP+n?=(rb), z2z=hC, (3.8) 


and eliminating z(z), v(z), w(z) and p(z) from 
(3.5)-(3.7), we have 


@ NS (oh? \? a bh ap 
-— 2p? 43 J 
ie : *) +( y ) Geum wares |r 


=0. (3,9) 


In order to get a quantitative result, we shall 
put 


(3.10) 


at €=0 and 1. (3.10) states that the tem- 
perature at the two free boundaries z=0 and 
h are kept constant. These boundary condi- 
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tions are satisfied by putting 
G= Iota neat. ) oll) 


Inserting (3.11) into (3.9). we get the follow- 
ing condition for the stationary convection 


4 FANVD oe 
ks a ian a ) 2, 
b v le 


Tosin sx€. 


(3.12) 


where 


4 
pet ORE 


re (3.13) 


In the case when the value of 2 in (3.13) is 
larger (smaller) than that given by the right- 
hand side of (3.12), the fluid motion is seen 
to be unstable (stable) and of developing 
(damping) character. If they are equal, that 
the 4 is a critical value which marks the 


stable and unstable motions. Keeping | 
Vv 


as a parametric constant and making use of 
(3.12), we can determine the critical value 2 
as a function of wave numbers 6 and s. The 
value of the critical 2 may become a minimum 
for certain values of 6 and s. As is easily 
understood, the mode of fluid motion cor- 
responding to these values of 6 and s is no- 
thing but the easiest one to be excited. Since 
A in.(3.12) monotonously increases with s, we 
have s=1 as the value of s for which 4 in 
(3.12) becomes a minimum. Thus putting 

di 
sal and db 
the value of b for which 4 becomes a minimum. 
The value of 5 thus obtained is an increasing 


ae 
y 


=0 in (3.12), we can determine 


function of e In any case, inserting 


the value of 6b thus obtained into (3.12), we 
get the minimum critical value of 4. This 
minimum value is once more an increasing 
function of ae Ms ‘ie 
existence of these relations is as follows. In 
convective fluid motions, temperature gradi- 
ent is an unstabilizing factor, while viscosity 
and Cortoutis force are stabilizing factors. 
The stabilizing action of Corionis force is 
investigated and well recognised in dynamical 
meteorology. In (3.12), the second term of 
the right-hand side of the equation denotes 


reason for the 


this stabilizing action of Cortouis force. The 
first term of the right-hand side denotes also 
the stabilizing action due to viscosity, while 
the left-hand side of (3.12) does the unstabiliz- 
ing action due to non-homogeneous heating. 
Equation (3.12) shows that the convective fluid 
motion can exist stationarily if these stabiliz- 
ing and unstabilizing actions balance with 
each other. Thus we can understand the 


2 (minimum)~(“°F") relation above _ ob- 
v 


tained. Furthermore, as is shown in dynamical 
meteorology, the stabilizing action of Coriouis 
force is more powerful for larger scale fluid 
motions, while the same action due to visco- 
sity is more powerful in smaller scale motions. 
These circumstances are reflected in (3.12). 
Thus, the first (second) term of the right- 
hand side of (3.12) becomes larger for larger 
(smaller) value of 6, that is, for smaller 
(arger) scale motion. In short, the action of 
Cortouis force makes fluid motions of smaller 
scales easier to be excited relative to larger 
scale motions. This is the reason for the 


existence of p( ee) relation above ob- 
y 


tained. From the above physical considera- 
tions, we may safely assume that the similar 


2 (min.)~(“8") and bea relations as 
v V 


above will be obtained under other boundary 
conditions than assumed above. The existence 


of the above eel relation is favour- 
vy 


able for our present study. The reason for 
this is as follows. 

Since S,2°-type fluid motion is rather of a 
smaller scale, we cannot expect that the S,2° 
motion becomes the easiest motion to be ex- 
cited (as is expected in (2) of §2) without the 


above pie relation. In short, while 
Y 


the S,2° motion is not the easiest one to be 
excited in the case when w=0, we may expect 


. 2h? - 
it will become so for certain value of ( ——}, 
vy / 


This is the most important result obtained in 
this section. 
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§ 4. In this section, we shall study about 
convective fluid motions under a uniform 
‘magnetic field. As is shown in Fig. 2, mag- 
netic vector of the uniform field is assumed 
to be in the yz plane, making an angle 0 
with the y axis. In this case, we have 


z 


Ho 
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curl H=4rJ, curl B= oe ; 
divH=0, J=o[E+uuxH)]. (4.2) 


In (4.1) and (4.2), # and o are magnetic per- 
meability and electric conductivity of the fluid 


respectively. u, H 2 E and i mean velocity, 
magnetic field, electric field and electric current 
vector respectively. From (4.2), we have 


OF cunt G (ux H) = d 
div H=0. (4.3) 
For H, we shall assume 
0 ie 
P| Eln “x0 ) ae a : (4.4) 
Fig. 9 Ay sin @ h, 
nm in which (hz, hy, hz) is a deviation of H from 
Bar =—prad pt+pypu—pyt+uJxH, the initial uniform field AH, and shall be 
treated in the same way as (uw, v, w), 0’, p 
Te, oT and J’ were in the last section. From (4.1)- 
Dt alae (4.1) . 
(4.4), we get the following equations. 
Ow aah Oh: Oh eo: ) 
—_ = )= — =. sitar i — p 
at ree s0( pili a \+ Lprw Peis ae etc., (4.5) 
Oh» Ou Ole 1 
—*¥—(0)= e __ p72 
oe Hy) cos er nat sin ss Pe Oey 5 TUCS (4.6): 
wie 2 
a 0=kp?T+ Bw , (a7) 
Ou , Ov dw _¢ Oh: , Oy  Ohe _ 
Ox Oy | Oz Ou ° Oy ° Oz Sa 
Eliminating wz, v, . from (4.5)-(4.8), we get 
T= o” oO ow 1 PET 9? ) : On\er 
pe gat eh eee a Pr. (4.9) 
Transforming (4.9) by the following equations 
4 
T= Tac es e=he m+n? = (xb) 
we have ie 
le —(xb)'] rand (xbP?T — ae oe e in cos 0+ af a 2 
ie sin ‘idl (San nov) T= (4.11) 
7 
Putting 0= 3 tentatively and assuming (see (3.11)) 
Tce sin sx€. (s= 1238. eae ). (4.12) 
we get the following condition for the stationary convection 
— a On 
es - (st bys 4p pS +8) , (4.13) 
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where = 
7 ag Bh! 
ky 
From this equation, we see that the existence 
of the magnetic field stabilizes the convection 
ow Zh 
Ov 


current and we get A(min.)~ and 


vel blitdid 
ov 


4 (min. ) ates y and p~(E relations 


y 
in the last section. 


§5. 


spherical vessel of radius a. 
rotating rectangular coordinates as is shown 
in Fig. 3, we have the following equations. 


z 


Fig. 3 
Du 2 anh 
grt —grad p+uy u+20(w Xu)—og , 
IDE 


ik (5.1) 
where u, @ and g denote the vectors of 
velocity of the fluid relative to the above 
coordinates, of rotation of the spherical vessel 
and of gravity respectively. A heat generat- 
ing source is supposed. The rate of heat 
generation is assumed such that, in the ab- 
sence of conduction and convection, the tem- 
perature would rise at a uniform rate 0. 
Thus in the state of no convection, we have 
_ OT) 


5.2) 
Oar ( 


=ky’?T)+0 5 


relations of the similar types as 


We shall now consider convective fluid 
motions which take place within a rotating 
Referring to the 


In a sphere of radius a, this gives 


O yee B 5 
ia gn? ij \—= Te ale : (5.3) 
where 
=e .8e) 
Crake ane Se 
is the temperature gradient at v=a. In (5.4), 


@ and ock are the rate of heat generation 
(per unit time and volume) and the heat con- 
ductivity respectively. In the same way as 
in $3, we get the following equations from 
(5.1)+(5.4): 


pu 1 grad P+20 xi)--agT=0, (5.5) 
»p_UrdTy__ 8B (5.6) 
(5.7) 


In (5.6), w#, means the radial velocity in the 
spherical coordinates (7, @, ¢) in Fig. 3. In 
these spherical coordinates, we have also 


u= (Uy, uo, Ud), o=(o cos /, —@ sin.0;.0); 
g=(9, 0, 0). (5.8) 


There are two types of velocity vectors u as 
follows which are known to satisfy (5.7). 
(H. Takeucai, 1950) 


Type I 
0 . 
ad Ox e 
Uy = Fr. m (7) i Waren +Gn ‘ m(7) y Wr yMy 
Uz - ie 
TU; = (nF nym rGn tn) Wn sm 9 
TU 0 
6 
00 
ale m Wn m ) 5.9) 
iene va as 
: sin # O¢/ 


where Wr.m is a solid spherical harmonic of 
of degree m and is expressed as follows 


Win Sn ms o)= r” Pr™(Ayem?, ( 10) 


In (5.10) Snjm is a surface spherical harmonic, 
and P,™ is a Legendre function. For the 
velocity vector of type I, we have also 


div u=0 
a n OM ree ats AC nm 5 14 3)Gnm | Wi m+ 
r ar dr 
(5.11) 
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a 0 
(VU) « - : 
(7°) y =fn,m(7) ti Waymt+9n,m(7) | Y Wn,m, 
(7°): oe P 
V1), = (fn m+P7'9n,m)Wnr,ms 
<> ) 
r(Puyo 50 
~ =Sam| 5 9 | Wom 5g) 
ae sin 8 O6 
in which 
aFy, DN ON En 
nym am aaa 2Gem 
J dr’ Pear 
@WGn m 2(n-+2) dG» m 
0.11 — : = * . eles) 
I 2 dr: Yr dr ( ) 


The other type of u satisfying (5.7) is given 
as follows. 


yoewil 
Ux Oy y Oz 
= OW n,m —_ OVWien 
Uy = Lin ih) qe Ties Zz Ox i ) 
U OW x 5m 0 Wn sm 
2 a 
Ox Oy 
Uy =0, 
1 OW 
— IE: 7 = 2,7 , 
ge : ee 0 Od 
OW n ™m 
= —LIn m : - alll ) 
ud m1) a0 (5.14 
div u=0. (5.15) 
(7’u),=0, 
i LE OVans 
4. =i “a - rm F 
(Veo ‘ Cs 0 db 
OW yr 
4 o— —T nm ) ai ) . 
(771) (r a6 (5.16) 
where 
Chan 20+) dL. 
Nh tes 1m . . 
pe ne age ae 


Taking (5.5) and (5.6) into consideration, we 
may assume 


Ihre irr a) Vises D= Tlie TW ane (5.18) 


for a of type I. Corresponding to wu of type 
II, we may take 

L=0, p=0. (5.19) 
The Laplacian of T in (5.18) is given as fol- 


lows 
VRS WV eine 
a ab Pm ees AT nm : (5.20) 
( dr r dr 
Accordingly equation (5.6) is satisfied by 
putting 


B 


tn > aes (Fn m+7?Gnm): (G2) 
ka 


Furthermore, if we omit the term wo xu) in 
(5.5), we can make equation (5.5) to be satis- 
fied by putting 


i Q(r" IT n,m) 


v(nfn mtr On ym) 


ty ings dr 
+ag7rTn,m=9, (5.22) 
v geste IIn,m=9, (6323) 
— 
Vioge Oh (5.24) 


The existence of the term 2(@ Xu) makes the 
present problem more complicated. A method 
of dealing with this complication will be 
considered in the section to follow. 

§6. For the sake of convenience, we shall 


denote x in (5.9) and (5.14) by tnm(Z) and 


Yinnh) respectively, and try to make (5.5) 
satisfied by putting 


u=S\ [ttn,m(D+tin,m(IDI. (6.1) 
The @¢ parts of all Un, m(1) and nel may 
be taken to be cos me@Gn=0. 1.20.2 oe 
The reason why we may assume this form 
may be understood by (5.5)-(5.20) in the last 
section. Next, the following transformations 
of variables are made 


== A 


i 
Fr nc )= - =F nl 7) 
a a” 


Gp mE) — Gn m1), 
Whe mn a= 1k. m(&), 


We shall denote also 


YB eS == Le m(X), 
Ie Has Lene). (6.2) 


Oy aE) 2n AF nm(E) 
a — =2Gawn 9) eige rs 
de te dé (€) 
DY fic ile eee With these notations, equa- 


tions (5.21)-(5.24) are rewritten as follows. 
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tom E)= — PO Ln Fy mE) + EG n(E (6.3) 
4 29 1 ae” Tn, i) iene 
Mn, m(&)+ Gn, m(E)|— 
[nfn,mE)+ E79 nm(E)I pe gE 5 & Lnm(E)=0, (6.4) 
i n(E)— lin. m(&)=0 (6.5) La ym(E)=0 F (6.6) 


Next, operating Un nD (which will be denoted by up(L) hereafter) to (5.5) and integrating, 
we get 


\; \; |." (6.5)-tis(r* sin 0drd0d$=0 . (6.7) 
0J0 JO 


: ed (6.1) into (5.5) in (6.7) and executing the integration with respect to 0 and ¢, we 
F have 


\ C(E)E"9(npF'p-+£°Ga)dé att D(E\EeF dE =0, 6.8) 
where 


coer [ 18} 43 Ena-%e171,(B,a) , 


DE)=(6.5)2-+~" 


© (ne(mne+1)10(8) ie rabna "8 Lala, 8) (6.9) 


In (6.9), the symbols (6.4)g and ee mean the left-hand sides of (6.4) and (6.5) with the 
suffix 8 respectively, and 


iiGi= " (Ps)? sin 046 ; (6.10) 
0 
(a,6)=\" Ps dP 8s 29d0 | (6.11) 
eT 
Tae B= The, ays ie Te PoPs +o oe sin Pecos 66 (6.12) 
P,=P™cos 0),  Fa=Frg,mp, --+-- (6.13) 


: “Similarly, operating up (I) to (5.5) and proceeding in the same way as above, we get 
1 

\ B(E\ERM e+ V Led E=0 , (6.14) 
0 


where 


E(E)=(6.6)p— soe [ne(e+1)L( BIE" "e [Fab(a, B)+ (neat Guha, B)}. (6.15) 


The symbol (6.6)g in (6.15) means the left- Equations (6.16)-(6.18) should be written down 
hand side of (6.6) with the suffix #, that is, for all possible #’s. In writting down these 
Ig. From (6.8) and (6.14), we get equations, we should keep the values of all 
C&H)=D(h=E(E)=0. (6.16) 2,’s equal to a fixed value a (see (6.1)). By 
eye : i -(6 btained, we can 
i t 5 th tions, we have b solving (6.16)-(6.18) thus o ; 
a anon : d determine F's, Ge, Ze, Ze and Ts. Comparing 


Hell), Moe” de(6:3 
ee ane) (6.16)-(6.18) to (6.3)-(6.6) and (6.17), we see 
ne dF p = pee -(mp-+3)Ga=0 , (6.17) that the complications brought out by Corro.is 
g de dé : Spmten... 4 20a" 
force are shown in the terms headed by 


b=" k pine p+ EGP) (6.18) “+a (6.9) and (6.15). In fact, if we put o=0 
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or all the (a, 8) and L(a, 8) to be equal to La, B) 


0, we have (6.3)-(6,6) and (6.17) in place of ay ba pane do 
(6.16)-(6.18). In view of these circumstances, = (a, B)+ne(ne+1) 5 Gea apee scl } 
we may call f(a, 8) and (a, 8) as coupling (6.19) 


integrals. Velocity vectors ws and we with 4nd the coupling integrals (a, 8) and 1a, B) 


(a, B), (a, B)<0 may be considered to be Vanish unless 1—mg=+1. The values of 
coupled with the action of Corionis force. J:(a@, 8) and L,(a@, 8) for several a and § are 


Now it can easily be shown that shown in Table I. 
Table I. 
Ne™ mp” if; In Ny” mp” I, Ts 
ce in ae ee : 
qo 90%) a a 20° 40 fe ‘ 
vo 
12 48 ; 5 
30-20 | 35 35 20 = 30 = = 
40 40 Ri 
30 40 pees 8 40 
63 21 40 30 7 ra 
50 8©= 40 4 el | 20 80 
99 33 40 60 tan s 
a 84 15 x 28 7 pre 
ov — . 0 La ete = 
: psy eee) eS a 11x13 11x13 
12 12 16x7 96 x7 
ie 21 = — 60 70 Se eS 
5 5 | 13x15 13x15 
96 384 
1 Ll } eee 
Aes 35 35 ee ees, wae 
200 200 | 4 12 
1 4. = == 21 | 
: 21 e 5 5 
320 640 / 192 384 
BL 1 — — OF 81 << 
tes We 33 11 ee 35 35 
135 x 4 735 x 20 40 200 
Bl 61 = ——- 1 a 
11x13 11x13 + ae) 7 7 
9 32 a 192 384 4l Bl = x 640 
7 as 11 
ee 480 480 x 5 6.5 ve x 20 735 x 20° 
7 7 TB IS8 ot | Tsai 
ae 480 x7 480 x 28 g om | _16x8x14 96x8x14 
4 4 ee Coe i 138 x5 13 x5 
ae 2800 x 24 2800 x 168 
2 pp Z I 
Lass sedis 96 
|| 256 x 63 256 x 63 x 6 320 Dp : = 
ae |  — a 7 7 
13 13 
: eure ae” _ 800 2400 
33 43 —800 15 x 160 7 7 
| "ig eee 82x70 480 x 28 
aan 82 x 82 x 35 32x 82x 385x6 11 11 
11 11 een 735 x 128 2800 x 168 
3 g«| «SC 8408024 315 x 40 x 1008 3 1x18 11x18 
11x13 11x18 ses 192 x 63 256 x 68 x 6 
7 63 120 x 72 x 56 960 x 72 x 56 ae 
iS 18 7 al Silex 64x 7 567 x 647 
7 gs | _220x126x72 220 72x 14x 63 i vi 
17 17 
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§7. As is stated in $2, we are most 
interested in the fluid motion of S,2°-type 
within the earth’s core. In the way of 
writing down in the last section, S.2°-motion 


is nothing but the Ue, oJ) motion in (6.1). On 
the other hand, it is shown at the end of the 


last section that the Uo, 2 (I ZT) motion cannot 


exist by itself under the influence of Coriouis 
force. It must be accompanied by the fluid 


motions of Us, AIT), tt1,2(1), Us,9( IT), ..types. 
These fluid motions are coupled with one 
another as is shown in (6.16)-(6.18). The 
way of coupling is shown schematically as 
follows 


ap ee Me ad) = uy oT) — us, x( 1) = PICO Ac 
(7.1) 


Similarly, we have groups of fluid motions 
as shown below 


1,01) 2 tz 0 LL) 2 u3,0(L) 2 tee (7.2) 
M1 01T) 2 te, 0(L) 2 3,0) 2 ..- (7.3) 
Mi iT) 2 te 17) 2 3,17) 2 --. (7.4) 
M101) 2 tayi(1) 2 gi) 2 ... (7.5) 
Uso IT) 2 3,211) 2 U2) 2... (7.6) 
Rem Vee UI) ) us g Iie «. < (7.7) 


ete: 


In general, there exist groups of fluid motions 
which contain w#,m(1 or IJ) with n=m as 
their first member. As is easily understood, 
the case when nm=m=0 is an exceptional one 
from this general rule. The problem (2) stated 
at the beginning in §2 is thus reduced to 
show that, under a certain condition, the fluid 
motion in (7.1) is the easiest type to be ex- 
cited among (7.1)-(7.7). 

Taking (6.16)-(6.18) and the results in Table 
I into consideration, we have the following 
results for the case in (7.1). In what follow, 
we shall omit the common suffix m=2, and 
Snell waite Simply’ iy, Gy, .... ‘for? Fe;3(€), 


Gy,2(&), At first, by (6.17) and (6.18), we 
have 

Zak, -dGs 

P aaeige 

Andhe dG, 

7G,=0, (b) 
£ dé + €& dé +7G, 
ied "SS" 8-4 a) eee (7.8) 


Glin alae pe ; 
— = PHAM “| a) 5 
dé + E dé (2F.+&G,) (a) 
Gilgen dT, _ bets j 
dé ie Pade k —=(40-2°G 4) (b) 
CICA o  —4) Pl acy gl ae ee (oS) 


Corresponding to C(E)=0 in (6.16), we have 


Caer + 2 ag PO) ae +e ae) 


1 dil.) +2 5 (- a 
—— ees T+ 
pvE dé Cat as os 
=0 ? (a), 
CFL. 8d; ) tae 12 ar il 
4 2G. 2( — 
(Ge te ae PO) +E Ge te ae 
il aE* Ts), AYoa 
er 1 
pe dé - 
480 4807 
—_-—— bea 10h ib), 
+7999 Le, as Ls) ©) 
Cte eee eee tr ee (7.10) 
rn C0) 
fipeplity (7.11) 
V 
and g is put as 
I= oF 5 (7.13) 
where gp is the gravity at E=1, that is, at 


the surface of the earth’s core. The gravity 
g in a uniform  self-gravitating sphere is 
shown to vary as is shown in (7.12). (The 
gravity in the earth’s core varies in a very 
similar way as in (7.12)). This is the reason 
why we have put y=%£. Now, corresponding - 
to D(£E)=0 in (6.16), we have 


da Fa 4 dF, TEs, (Dany O88 oe 
2 Ls NG) 27,,=0, 
dE ne dé +2G,— Fai ae 7 &*L; 
(a) 
ae Fy, 8 dF, Ih 
a 2G 
ge de 
1 480 x5 AbAQ) esi )- 
Th 2 =—=(0). b), 
+5552 7 3+ it E (b) 
ics OP VN NS eo 8B (@/oiles}) 


Lastly, corresponding to E(€)=0 in (6.6), we 
have 


PL, 8 aly _ 7 afer, 
de t & dé 1440 7 

we 
+ Moreen [+ R44, 


7 
+260 | =0, (a) 
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lig, edly aM g pnd OAT, 
de? ' — d— 210x120 


2800 x 24 
pple Ng +660) |+] aoe 


11 11x13 
2800 x 168 Soa 
SOE GEG |e 0 b 
apie >| Y 
gt ge eae ee (7.14) 


If we omit wp,. (J or I) with n=4_in (7.1), 
we get the results (a) in (7.8)-(7.14). Thereby, 
the F, and G, terms in (7.14) (a) must be 
omitted. The results thus obtained will be 
called the (a)-approximation. Similarly, the 
results obtained by taking the (a) and (b) 
equations in (7.8)-(7.14) and thereby omitting 
the Fs and Gz terms in (7.14) (b) will be called 
the (b)-approximation. In the (a)-approxima- 
tion, we must solve five differential equations 
among five unknown functions (F., Gz, I], To, 
I;). In the (b)-approximation, the number 
of unknown functions increases to 10. In any 
case, proceeding in this way to (c), (d), ...- 
—approximations, we shall obtain the required 
solution for the case in (7.1). Having thus 
established the fundamental equations, our 
next thing to do is to consider about the 
boundary conditions in this case. By the 
assumption of ‘‘ viscous ’’ fluid, we must have 


Unym(L), UnmD)=0, at E=1. (7.15) 
As the boundary thermal condition in this 
case, we shall take rather tentatively 

ili) aL een (eS) 
By (7.16) is meant that the temperature at 
y=a is kept constant. This is a plausible 
condition to be taken for the earth’s core 


boundary. Taking (5.9), (5.14) and (5.18) into 
consideration, we have 


ME ri nt? Gn ,m=Gnim=Ln,m= Tn, m=0 
clare CUNEO 
Frnim, Gnim, ---. must also be finite at r=0 
(or €=0). Thus our problem is reduced to 


solve (7.8)-(7.10), (7.13) and (7.14) under the 
boundary conditions in (7.17). It is easily 
seen that our problem thus established is an 


: v : 
eigen-value problem for 2= © Boa contain- 


v 


a 


ing a parameter Q=°° , Although the above 
VY 


results are obtained for the case in (7.1), simi- 
lar results may be expected to be obtained for 
(7.2)-(7.7) equally. It is, however, almost 
hopeless to solve the eigen-value problem thus 
obtained by trial and error methods. In order 
to solve the problem in its simplest form, 1.e., 
in the form of the (a)-approximation cited 
above, we must solve five simultaneous dif- 
ferential equations of second orders. A 
method of dealing with these difficulties will 
be devised in the next section. 

§8. Taking (7.17) into consideration, we 
shall put 

NE n m+ &?Gn,m= (1 —E&)(Ay +A,E+ O48 as 

Lnm=(1—E)( Bot Bié + tee 2 
where Ay, By, Ai, Bi, ..-. are undetermined 
constants. In order to avoid unnecessary 
confusions, we shall denote Ag(2, 22), Bo(m, 2), 
.... for the above Ay, By, .... The equation 
(6.17) is transformed into 
@(nF yn m+E?Gn,m) 
Ed= 

The value of Gn,m is calculated by (8.1) and 
(8.2). Inserting Gr,m thus obtained into (8.1), 
we get Fuim. The form of 2Fnj,m+&?Gn,m in 
(8.1) is so adjusted as to make Fa,m thus 
obtained vanish at €=1. Next, inserting (8.1) 
into (6.18), we have the following differential 
equation for Tym. 

@Tam | 20+1) dTrm__ Ba® 

d= £ GER NE 

where A;’ is a certain constant which is a 
function of Ao, Ai, .... in (8.1). The equation 
(8.3) is solved in the form 


Pan) = Tees ee 
OR “TF (j+2\ G+3+42n) 


(8.1) 


(72 + 1G, ead 


(8.2) 


BAYE, (8.3) 


EI? 
(8.4) 
In (8.4), A’ is a constant, which is deter- 


mined so as to make Ty,» in (8.4) satisfy the 
condition (7.17). The value of A’ thus ad- 


justed is 
= Ay 
GLDG+S4 In) ee) 
Inserting (8.1) and (8.2) into D(E)=0 in (6.16), 
we can calculate In». Thus, with 


MP nmt+E Grim and Lnim in (8.1), we can make 


ont A 
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three groups of differential equations (6.17), 
(6.18) and D(E)=0 satisfied among the required 
five. There remain, still, two groups of dif- 
ferential equations to be solved. They are 
Cé)=E(—)=0 in (6.16). We have also two 
groups of undetermined constants Ay, A, By, 
B,, .... We shall now try make C(é)= F(6E)=0 
satisfied approximately by choosing Aj, Aj, 
-.+-; Bo, Bi, .... adequately. Inserting (8.1)- 
(8.5) into (6.14) and 


se 
[ CHEM A e+ EG dE=0 


(see (6.8) and the equation D(é)=0 above 
satisfied), we have equations of the following 
forms 


1 
As C(EE"0(1—£)°dE=0 
0 


(8.6) 


A,\ CEE EdE=0 ; (8.7) 
B,|) BEC —BdE=0 . 
Bl BOE Seee 0. (8.8) 


These equations will be used to determine 
the values of Ay, Ai, By, Bu, In order 
to show the way of determining Aj, A, .... 
by an example, we shall take up the (a)-a 
proximation for the case in (7.1). C(&) and 
E(&) for this case are shown in the left-hand 
sides of (7.10) (a) and (7.14) (a). (In (7.14) (a), 
~ F, and G,; have to be put equal to zero). We 
shall, furthermore, take only the first terms 
of 2F,,.+&G,,, and Z,,, in (8.1). The results 
obtained in this way will be called the (a, 1)- 
approximation. Similarly, the results obtained 
by taking the first two terms of 2F,,.+&?G»,» 
and Z,,. will be called the (a, 2)-approximation. 
fee Wy a(Ds) 2),,.AC ek )etc, 2), . approxi- 
mations may also be defined similarly. The 
equations to determine A, and B, in the (a, 
1)-approximation are obtained by (8.7) and 
(8.8) as follows. 


Sai we 2A 
ii 18-50-.21e28e15 5) 
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In order that these equations are compatible, 
the determinant formed by the coefficients of 


A, and By must be equal to zero. The 
equation for 2= cad thus obtained is 
A= 10602.5-+1.96994.0, 9-9 19) 
7) 


It is to be noted that this is of the same 
form as in (3.12). The values of 2 in the 
(a, 1)-approximation are calculated for several 
#° and are shown in Table II. In Table II, 
are shown also the values of 4 in the (a, 2)- 
approximation for the same 2?’s. The values 
of 2 in the approximation differ but little 
from those in the corresponding (a, 1)-approxi- 
mation. The (a, 3), (a, 4), .... approximations 
are expected to give results practically the 
same as those in Table II. 


Table II. 
Q ie Canty (a. 2) 
0 | 1.060-10! 1.06-104 
102 | 1.079+10! 1.08-10# 
103 1.247-104 | 127-104 
10! | 2..920-10# | 303-10! 


Thus, in order to get the value of 4 in the 
lim (a, j) approximation, we have only to 


jroo 
follow the (a, 1)-approximation. Similarly, 
Pie lime (by. isa ith Cera), . approximations 


Paes 
of 2 may safely be inferred from the (b, 1), 
(Gee 1). . approximations. The values of & 
in the (b, 1) and (c, 1)-approximations for the 
case in (7.1) are shown in Table III. In 
Table III, are shown also the values of 4 in 
the (a, 1)-approximation for this case. The 
(a, 1) value is taken from Table Il. The (c, 
1) value of 4 in Table HI differs but little 
from the corresponding (b, 1) value. The (d, 
1), (e, 1), ...-approximations are expected to 
give results practically the same as those in 
Table III. Thus, we may use the (c, 1) value 
in Table III for the exact A in the case (7.1). 
Similarly, the eigen-values for the cases (7.2)- 
(7.7) are calculated and are shown in Table 
IV. A discussion on the results in Table HI 
and IV will be made in the next section. 
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Table Il. i (min. ~(79"") relations in §3. From our 
——— v 
2 (a, 1) (b, 1) fe, A) present point of view, it is important that 
at, Fi 060 x10! | | we have the quasi-S,2° motion which becomes 
102 arent | 1.078104 | 1.078x10! the easiest type to be excited in the case 
108 (ise ene 1.23610! | 1.235x1lot When 2?--10*. Since our previous studies in 
104 See. 2.259x10' | 2.311x10t Paper I, If and II were made under the 
assumption on the existence of S/° motion, 
Table IV the result obtained above is very favourable 
2 ee to us. Though the existence of T;?, S.2, ... 
2 uo) u1,oUT ) uyL) motions in the quasi S-motion may bring a 
MMMM Oo) 1s § us ek ca oe ae , little alteration to our results obtained in 
jt be 0.811 x 10 : 
. ona! re fee Paper I~III, the alteration will not be so 
+ ee che acre large. We shall now estimate the value of 
Me oe ie kinematical viscosity » corresponding to 
104 2.6 3 5 ~ 
he See & be 
— : aoe 0% Q= ==10? 9.1) 
wD) us,a(ID) uta(Z) ins aa . ( 
0 1.06x10 | 1.5810 | 1.58x10 Putting o= se (secyok and naes3 4 Over 
102 1.10 1.55 1.54 : 
10° 1.48 1.779 1.64 into (9.1), we have v-=10"5" _ In short. if 
sec 
10¢ SAr 3.562 2.44 ‘ ; ‘ oe 
the kinematical viscosity of the fluid in the 
§9. In this section, we shall call the mo- garth’s core is about 10° , so we have 


HOMSMITING/ el yn /e2)) (1c); eee ras Aut MmoOcons 
of quasi-S,2°, quasi-S;°, quasi-S,°, . types 
respectively. The teason why we call in this 
way is as follows. As was stated at the 
beginning of §7, the w.,,(7) motion is nothing 
but the S.?° motion in the sense of Paper I, 
II and III. It is usually called the poloidal 
motion of m=m=2 type. Similarly, the 
m#,0\1) motion is the poloidal motion of 2=1, 
m=(0 type. It is usually denoted as S,°. As 
is seen from Table III and IV, if we put 
12MIN (ls (ied) putea) sake, Wier Dave. the 
eigen-values 2=1.06 x10‘, 0.811 x 10‘, 1.06 x 104 
The eigen-functions corresponding to 
these eigen-values are of pure S,2°, S,°, S,°,. 
... types. (H. Jerrreys and E. M. Buanp, 
1951: S. CHANDRASEKHAR, 1952). In any case, 
we shall now study about the problem stated 
in (2) in §2. The minimum values of 2 for 
the cases 2?=0, 107, 10° and 10‘ are given 
by the fluid motions of the quasi-S,°, S;, S,}, 
and S,? types respectively. The minimum 


A becomes larger for larger 2?. These rela- 


2 2 
tions correspond to the p= (= and 
Vy 


the fluid motion of (quasi) S,?°-type there. 
The viscosity v in the earth’s core was hitherto 


estimated to be v-==10-2~10-3—™_ (Ee Ge 
sec 


BuLLaRD, 1949) This seems contradictory 
with the result above obtained. It is the 
opinion of the present authors, however, that 
this contradiction is not so fatal to our 
studies. What our present study is concerned 
with is the eddy viscosity, whereas the vis- 
cosity y=10-?~10-* above referred to is the 
molecular one. It is a well-known fact in 
dynamical meteorology and oceanography 
that the eddy viscosity for large-scale motions 
is far greater than the molecular viscosity. 
We shall now think over the values of 2 ob- 
tained in the last section. The lowest value 
of 4= aabige in the case 2?==10' is that for 
the quasi-S,2° motion. 
is shown to be 


In Table III, this value 


A-=2.3x10!. (9.2) 


Inserting g=10°, a=3.4x108, ~=10-5 and 
v=10" into (9.2), we have 
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Doe i0- 


; (9.3) 


For fluid motions in the earth’s core, the 
thermal diffusivity Rk in (9.3) should be inter- 
preted as the eddy diffusivity. In contrast 
to this, the thermal diffusivity Rk in (5.1)- 
(5.4), or in the state of no convection, should 
be interpreted as the molecular diffusivity. 
The molecular & in the earth’s core is esti- 
cm? 
sec — 
lecular F is a little greater than the molecular 
v, the eddy k& is considered to be a little 
greater than the eddy »(=10'!). Thus inserting 
k=10"~10” into (9.3), we have B=10-5~10-°. 
As was stated before, 8 is the temperature 
gradient (at the surface of the earth’s core) 
in excess of the adiabatic one. The latter is 
estimated to be 10-°. Thus it is shown that 
the convection current of quasi S,?°-type can 
exist stationarily under the thermal gradient 
which is of the order of magnitude of the 
adiabatic one. This is the thermal state 
which is usually inferred to exist within the 
earth’s core. Thus the results obtained in 
this section may be summed up as follows. 
In order that the covection current of the 
quasi S,?°-type exists stationarily in the 
earth’s core, the kinematical viscosity v there 
must be of the order of magnitude of 10". 
“If »y is taken to be 10", then the thermal 
diffusivity k there must also be taken to be 
1011~10". If these values for v and k are 
admitted, the thermal gradient there is shown 
to be of the order of magnitude of the adiabatic 
one. This is the thermal state which is 
usually supposed to exist in the earth’s core. 
§10. Since the value of v obtained in this 
paper is rather too much greater than is sup- 
posed hitherto, details of the tentative con- 
vection current model in §5 of Paper II will 
have to be changed a little. The magnetic 
field and the convection current velocity in 
the earth’s core must be taken to be a little 
smaller than supposed there, but the altera- 
tions will not be so large. In the discussions 
in §5~-9 of the present paper, the existence 
of the magnetic field was not taken into ac- 


mated to be k=10- Since the mo- 


count. The existence of the magnetic field 
within the earth’s core, however, may bring 
no drastic changes to the results obtained in 
the present paper. This is seen from the 
equation (4.13). The result in (4.13) was ob- 
tained for plane boundaries and those for 
the boundary conditions other than those in 
§ 5-9, 

By the same reasoning as that at the end 
of §3, we may safely conclude that the quasi 
S.°°-motion becomes the easiest to be excited 
within the earth’s core in the case when the 
second term in the right-hand side of (4.13) 
becomes approximately equal to the first. In 
calculating the first and second terms, we 
may take also s=b=1. Thus we get the 
following condition for the existence of the 
S.?"-motion with in the earth’s core 

ow Heh? 
pv an 
Pemoine GaN mail, weaSA ate, o= ile 
and v=10" in (10.1), we have Aj==10. Since 
the result in (4.17) is obtained by putting 


50 (10.1) 


7 


t= 9 in (4.4), the value of Hj above obtained 


is considered to be related to the radial com- 
ponent of the magnetic field within the 
earth’s core. The value of the above Mh 
interpretted in this way agrees well with that 
which has usually been inferred for the earth’s 
core. 
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Abstract 


Wave phenomena which are seen along the surface of a semi-infinite elastic solid 
when subjected to an impulsive line force in the direction normal to its surface have 


been studied. 


The results obtained are as follows. 


(1) Immediately after the impulsive force is removed, the vertical displacement 
is downward in the neighbourhood of the origin, and it is upward outside that down- 


ward domain. (§3) 


(2) In the downward (upward) domain, the horizontal displacement is toward 


(away from) the origin. (§ 4) 


(3) The surface deformations thus produced are propagated with the velocity of 


RAYLEIGH wave without conspicuous changes in form. 


($3 and § 4) 


(4) The wave length of the RAYLEIGH wave is proportional to the time of duration 
of the impulse, and is of the same order of extent as that of the downward domain 


at the moment when the force is removed. 


§1. Wave phenomena which are seen along 
the surface of a semi-infinite elastic solid 
when subjected to an external impulsive force 
were studied by H. Lamps (1904). His results 
are chiefly concerned with the phenomena 
“which are to be observed at the places far 
distant from the wave origin. The deforma- 
tion of the elastic solid in the neighbourhood 
of the origin was not studied in his paper. 
On the other hand, an experimental approach 
to this problem has recently been made by 
K. Kasawara (1952, 1953). In his experiments, 
an apparently curious phenomenon has been 
found in the neighbourhood of. the origin. 
He has stated: “When an impulsive vertical 
force acts downwards on the surface of an 
elastic medium, a part of the surface within 
a certain region begins in upward motion. 
The disturbance which occurs in such a way 
is propagated as time passes on. Such a 
phenomenon seems rather curious from the 
stand point of the theory of the statical de- 
formation of elastic bodies. Because, as 


($ 3) 


Bousstnesq has studied, every part of the 
surface of a semi-infinite elastic body sub- 
jected to an external vertical force acting 
downwards at a point, is displaced downwards 
and no part showing upward displacements 
can exist on the surface. In fact, a statical 
experiment carried out on the same medium 
produces such results”. KAsaHARA has shown 
theoretically the possibility of swelling up of 
the surface of an elastic solid when subjected 
to an impulsive downward force, but has 
failed to bring a unified view in connecting 
this upward motion with the phenomena ob- 
served at the places far distant from the 
origin (1953). It is the purpose of the present 
paper to show a simple way of calculating 
the surface displacement of the elastic solid 
at any time and at any point and to fill up the 
gap between the knowledges concerning two 
extremes above cited. 

§2. It was shown by Lams (1904, 
equation (93)) that the vertically downward 
displacement v at (#, y=0, z) due to the 


his 
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vertically downward force Q(¢) at the origin (w=0, y=0) is given as follows. 


a NA 
pq (262—b?)' +1600? -—@)\(b°— 
eile. \ me —x0t—02)00 (2.1) 
zu” \y (26°—b?—402\/ 0 — Sai 


In (2.1), is the distance of the point from 
the origin and ¢ is the time, @ and 6 in (2.1) 
are 


De Be at les, 
Oat ' Ve 
respectively and the symbol p denotes the 
principal part of the integral which follows 
it. It is to be noted that the equation (2.1) 
is valid for any x and ¢. Equation (2.1) can 

be rewritten as follows. 


wtrennaf (lh awa) 
Ei \S (a é) a{ Tate) 
(2.3) 


(2.2) 


in which V(@’) =V(@/a) is 


a m3(20’2—m?)\/ 672-1 

Se leg (20°7—m?)' + 160"(0—1 (9? — 02) 
for 1<0'<m 

vale m/ 0/2—1 

Vo )= (20°°—m 249”: /69/2— ly ie 29992 

for ma 
tte: & Ve (2.4) 
a Veg 


The value of V(@’) for the case 4=y is given 
in the Lams’s paper. We reproduce it in 
Laem Ae 

We shall now calculate the value of wv) in 
the case when 


Q(t)=const@, for —r<t)<c 
and Q(t))=0, otherwise. (2.5) 
Putting f=ct tg = Cty (2.6) 
in (2.4), we have 
Vere 41) 
5 EEE p| i (—V(0"))dé’. (2.7) 
q@ YER yh ==i1l)) 


This is the formula by which the value of 


Bis. 1. 


The function V(@’) 


Vo at any x and ¢ can be calculated. The 
hatched area which is shown in Fig. 1 re- 
presents the integral on the right hand side 
of (2.7) corresponding to a certain 2’. In 
executing this integration for various values 


of z#’, we come across a difficulty when the 
domain of the integral Yea includes 
6’ =1.88.... — At @=1.88... the integrand 


—V(9’) in (2.7) becomes infinite of the order 
of (6’—1.88...)-1. This difficulty is similar 
to that encounterd in executing the integral 


+1 da 
> 
a) 


by the aid of the principal value expressed 
the symbol # in (2.7) removes the difficulty 
stated above. Thus making use of (2.7), we 
can calculate the value of —v) (upward dis- 
placement) for any x and ¢. The results of 
calculation for the case when A=. are shown 
in Fig. 2. 

§3. As is seen in Fig. 2, the vertical dis- 
placement at ¢=0 or t=t is downward for 


for instance. The limiting process 


v x : : 
av <0.49 or V <0.98 respectively. It is of 
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Figs 2. 


interest to note that the linear extension of 
the downward domain is proportional to rt 
(half the time of duration of the downward 
force). Thus if we make t— oo, we have the 
downward displacement everywhere on the 
surface of the elastic solid. This is nothing 
but the result obtained by Boussinese in the 
case of a statical downward force. There is 
a domain of the upward displacement outside 
the downward domain above stated. In the 
two dimensional problem as ours, the points 
of zero displacements lie on a straight line 
x=const. at a fixed time. As time passes on, 
this line of zero displacement, i.e., the 
boundary line between the upward and the 
downward domains is propagated from the 
origin with the velocity of RAYLEIGH wave. 
On the other hand, the front of the upward 
displacement is propagated from the origin 
with the velocity of P wave (compressional 
wave). The origin times of these two “waves” 
are both ¢=—r. Thus the upward (and also 
of the downward) domain increases in its 
extent according to time. The peak of the 
upward displacement is propagated from the 
origin with the velocity of RAYLEIGH wave. 
The origin time of this upward peak is also 


5 6 


( Distance Parameter; 


non Dimension ) 


Upward deformation —vy 


t= —t, i.e., at the time when the vertical force 
begins to act. The trough of the downward 
displacement is propagated from the origin 
also with the velocity of RAYLEIGH wave. 
The origin time of this trough is #=+7, i.e., 
at the time when the vertical force is re- 
moved. The time interval at any point be- 
tween the arrival times of the upward peak 
and the downward trough is 2c. The corre- 
sponding “wave length” Z of the RayLEicH 
wave is Vr-4r, that is, the RAYLEIGH wave 
velocity times twice the time of duration of 
the downward force. As was stated at the 
beginning of the present section, the domain 
of the downward displacement at ¢=T is 
characterised by the relation |#|<0.98cVe. 
and the width D of the domain of the down- 
ward displacement is 2x0.98cVp. Making 
use of the relation Ve=0.53V>, for the case 
when A=p, we have D=3:7tVr. This is 
approximately equal to the wave length 
L=4tVpr of the generated RAYLEIGH wave. 
In Fig. 2 are shown only the displacements 
up to the time ¢=5r. The approximate pattern 
of the displacement for ¢/>5, however, can 
be obtained as follows. The front of the up- 
ward displacement is at x=V,r(t+r). The 
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peak of the upward displacement and the 
trough of the downward displacement are at 
w=V,(t+rt) and Vx(t—r) respectively. The 
displacement in the neighbourhood of the 
origin decreases gradually as time goes on. 
In other respect, the pattern of the displace- 
ment is similar to that at fdr. 

The curve of displacement at z=0 is shown 
in Big. 3. 


Figy 3. 


Displacement at the origin 


For the large values of t=r?’, the curve is 
expressed approximately by the equation 
V 374 Sys 2 
_ V)=——— log, 1 : 
gq wae toe(14+ 7) 
The equation (3.1) is obtained as follows. 
From (2.4), we have 


(3.1) 


—ViO)= (3.2) 


ie 
4 0 
for the large #’. Inserting (3.2) into (2.7), we 
get the equation (3.1). It can readily be seen 
that all of the numerical results above ob- 
tained are closely connected to the infinity 
point at #’=1.88... of the function —V(6’) 
in Fig. 1. This infinity point at 0’=1.88... 
is, in turn, closely connected to the RAYLEIGH 
wave velocity for the case when 4=y. In 
fact, the numerical factor 0.53 in the relation 
Vr=0.53V,~ is nothing but the reciprocal of 
this 1.88.... It is therefore expected that the 
similar results as those obtained above will 


be obtained for the elastic solid for which 
j-<y. KASAHARA’S experiment referred to at 
the beginning of the present paper is made 
with an elastic solid (agar-agar) for which 
4=4~5m. Adopting this value of 4/z, his 
experimental! results agree quite well with the 
above theoretical results. According to Kasa- 
HARA, for example, the peak of the upward 
deformation in his experiment is propagated 
with the velocity of about 1cm/5m sec 
=200 cm/sec. (see Fig. 4 in his paper IT) and 
this is seen from Table II in the same paper 
of his to be about the velocity of the RAYLEIGH 
wave in his experiment, just as is expected 
from our calculations. 

§4. By the similar methods as adopted in 
the preceding sections, we can calculate the 
horizontal displacement (2, y=0,7). The 
equation by which the value of mw is to be 
calculated is as follows 


27 f P 4 
"Qt —maHa(t Si = ) 


Vpt 
Veto 41) 
+2 a UO) de’, (4.1) 
Yeter—1) 


in which 0 is the Dirac’s function, and 
BV (20° —n?)V 077-1 V m0? 


Fg ema 6 i 
UO) = — “96735 — mi +1600? —1\ m0) ° 
1<0’<m 
u(0’)=0, 6’<1, 6'>m 


H=0.12500 (A=,). (4.2) 


The equation (4.1) and (4.2) are derived from 
the equation (92) in the Lamp’s paper, and 
the value of U(@’) for the case when 4=yz 
was calculated also in the Lamp’s paper. Thus 
making use of (4.1), we can calculate the 
value of mw) for any z and z. The results of 
calculation for the case when 4=y are shown 
in Fig. 4, 

As is seen from Fig. 4, the horizontal dis- 
placement at ¢=0 or ¢=t is towards the origin 
for x less than about 0.57rVp or 1.15rV> re- 
spectively. This domain of inward displace- 
ment roughly corresponds to the downward 
domain referred to at the beginning of the last 
section. Outside this domain, there is also a 
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| 
(Distance Parameter; 
non Dimension ) 


0.8 
Fis. 4. 


domain in which the horizontal displacement 
is outward from the origin. In the “ outward 
domain”, the vertical deformation is upward. 
The front of this horizontal displacement is 
propagated with the velocity of P wave. The 
origin time of this wave is f=—r. The 
principal part of the horizontal displacement 
at ¢ is limited in the interval 


Ve(t—t)<a2<VeR(t+r). 
“The value of ww in this part is about 


—0.7 nf 29 The front and rear of the 


principal part are propagated from the origin 
with the velocity of RAayLEIGH wave. The 
origin time of these front and rear waves are 
t=—rt and ¢=rT respectively. 

§5. Throughout the present study, a line 
source of wave has been assumed, and the 
problem has been treated as a two-dimensional 


Horizontal Displacement 2 


one. In the two-dimensional problem as ours, 
it is to be noted that the waves are propa- 
gated without conspicuous changes in form. 
In a three-dimensional problem with a point 
source of force, this is not the case. This 
latter problem will be studied in another 


paper. 
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Determination of Elastic Wave Velocities in Rocks by Means of 


Ultrasonic Impulse Transmission. 


By 
Akira Kusorera. 
Geophysical Institute, Faculty of Science, Kyoto University, Kyoto. 


Abstract 


Velocities of longitudinal and transversal elastie waves in rock specimens haye 
been determined from the transmission times of ultrasonic impulse (50 KC—800 KC) 


through them at ordinary temperatures and pressures. 


POISSON’s ratios have also been 


caleulated from the velocities of the two kinds of elastic waves. 
18 specimens of granite and 34 of marble, all collected from various localities of 


Japan, were examined. 


In addition, 12 apecimens of marble from Italy and various 
other kinds of rock, such as andesite, sandstone, tuff, ete. were also examined. 


The 


experimental results are summarised as follows: 


longitudinal | transversal POISSON’S 
wave (km/sec.)} wave (km/sec.) ratio 
Granite | 4.09-5.89 2.51-3.55 ; 0.19-0.28 a 
Marble 8.75-6.94 2.02-3.86 0.18-0.35 


Within the frequency range of our experiments, dispersion phenomena have been 
found in none of the rock specimens examined. 


Introduction 


§1. 

The elastic constants of rocks for low and 
audio-frequencies have been measured by a 
number of workers. In experimental lobora- 
tories, the constants can be determined from 


the resonance frequency of a rod of rock 


specimen which is subjected to standing 
elastic vibrations (K. Ima: 1939). In fields, 
on the other hand, the constants can most 


directly be determined from the time of 


transmission of artificial or natural seismic 
waves. If one wishes to use this direct 
method in laboratories, an inconveniently long 
specimen of rock will be needed, because the 
specimen must be at least several times as 
long as the wave-length of the elastic wave 
used. This difficulty will disappear if ultra- 
sonic frequency is used instead of audio- 
frequency. In the present paper, the writer 
is going to describe the results of his lobora- 
tory determinations of elastic constants of 
rocks from the time of transmission of ultra- 


sonic waves (50 KC-300 KC) in rock specimens. 
The dimensions of the specimens are 20cm x 
17cmx6cm in case of granite and 12cmx9 
cm xX(2~1)cm in case of marble. 


§2. Apparatus and Experimental Method 


Fig. 1 shows the general view of the ap- 
paratus used and Fig. 2 is its schematic di- 
agram. The apparatus consists of a pulse 
generator, pulse amplifier, sweep circuit, time 
mark circuit and cathod-ray oscillograph. 

BaTiO; crystals were used to generate and 
receive ultrasonic waves of 50KC, while 
quartz crystals were used for 150 KC and 300 
KC. For the latter case, X-cut crystals were 
used for longitudinal waves, while Y-cut 
crystals were used for transversal waves. D. 
S. Hucues and J. H. Cross (1951) used X-cut 
crystals for both kinds of waves, but accord- 
ing to the writer’s experiences, Y-cut crystals 
are more adapted for transversal waves. 

The pulse generator (including power am- 
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plifier) which excites the crystals has been so 
designed that the pulse is sent into the speci- 


men a short time after the sweep circuit- 


bigins to operate. Sweep signals are fed to 
the horizontal plate circuit of the cathod-ray 
oscillograph while the imput impulse to the 
rock specimen and the output impulse from 
it are both fed to the vertical plate circuit of 


Fig. 1. General view of the experimental 


apparatus. 


Oriving 
Crystal 


Rock 
Specimen 


« feceiving 
Crystal 


Cathod Ray 
Ooct)lograph 


Fig. 2. 
experimental apparatus. 


Schematic diagram of the 


the oscillograph together with time marks. 
The time marks are given at every 3-40.06 
microseconds. The apparatus was designed 
by Assistant Professor F. Konno of the Depart- 


ment of Electrical Engineering, Kyoto 


University, and was constructed by his as- 
sistant, Mr. T.. YAMADA. 
§ 3. Rock Specimens 

The rock specimens examined are as 


follows: 


Granite 18 Japan 
Marble Bo Japan 
Marble 12 Italy 

Marble 1 Korea, 


all collected by Mr. S. Nacaoxa. These were 
placed at the writer’s disposal by the kindness 
of Professor S. Marsusnira. The localities. 
from which they were collected are listed in 
Tables I and II. Each of the rock specimens 
has its own characteristic appearance and has 
the proper name according to-its place of 
occurrence as also listed in the tables. The 
properties of these specimens of granite and 
marble are believed to be such that they well 
represent those of the corresponding rocks 
which occur in their respective places. Un- 
fortunately, however, chemical and mine- 
ralogical compositions of the specimens have 
not yet been determined. 


$4. Results 


The cathod-ray oscillograms due to long- 
itudinal and transversal waves in rock speci- 
mens were taken photographically. Fig. 3 is 
one of the records obtained for a granite 
specimen. On the first record are shown the 
incident and received wave pulses due to the 
X-cut quartz crystal oscillations, while the 
third those due to the Y-cut crystal oscilla- 
tions. The second and fourth records are 
time marks corresponding to the records above 
them. As is clear from the records, when 
the X-cut crystal is used, only one phase can 
be seen in the wave trains received. This 
phase corresponds to the longitudinal wave. 
In contrast to this, when the Y-cut crystal is 
used, the wave trains received show two 
sharply distinguishable phases, the first 
smaller one corresponding to the longitudinal 
wave and the latter larger one to the trans- 
versal wave. It is likely that in case of the 
Y-cut, while the main oscillations of the 
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(1) Longitudinal wave of X-eut erystal. 

(2) Time Mark. 

(3) Longitudinal and transversal waves of 
Y-cut erystal. 

(4) Time Mark. 


Fig. 3. An example of records obtained for a 


granite specimen. 


GRANITE 


mers eer Sa 607 
Up AM%ee 


Fig. 4. Longitudinal and transversal wave 
velocities in granite specimens. 
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Longitudinal and transversal wave 
velocities in marble specimens. 


| hea 


crystal are in X-direction and produce the 
transversal wave, those in Y-direction due to 
the thickness oscillation produce the long- 
itudinal wave. 

The velocities of longitudinal and trans- 
versal waves obtained for various specimens 
are listed in Tables I and II for granite and 
marble specimens respectively. Figs. 4 and 
5 have been obtained by taking the transversal 
wave velocity vs in ordinate and the longi- 
tudinal wave velocity vp in abscissa. The 
points in these figures are seen to lie ap- 
proximately ona straight line within a certain 
range of velocity values. The velocities otf 
elastic waves in granite and marble range as 
follows: 


Vv» (km/sec.) | vs (km/sec.) 
Granite | 4.09-5.89 | 2.51-3.55 
Marble 3.75-6.94 2.02-3.86 


ee ee 
Porsson’s ratios as well as densities of the 
rock specimens are also given in Tables I and 
II. The range and the average values of 
Porsson’s ratio are as follows: 


| range | average 
Granite 0.19-0.28 0.23 
Marble 0.18-0.35 Ona 


ee  EEEEEEEEEEEE 
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Table I. Locality, Proper Name, Longitudinal and Transversal Wave Velocities, 
Porsson’s Ratio and Density of Granite. 
Pr | POISSON’S Density 
No. Locality 3 Rh Vy km/sec. | Vs km/sec. | ratio er /ems 
ie oes Mis ponies oa } - 

1 | Hyogo-Pref. | Akao-ishi | tee 3.27 26 2.86 
2 ” | Hon-mikage (1) 5.77 3.50 Aad | Pat 
3 ” ” ah 5.38 3.16 | .28 | ats) 
4 | Osaka-Pref. | Kitazima-ishi 5.89 3.37 | .26 | 2.76 
5 ” a | 5.63 3.16 Ben, Zieslal 
6 Aichi-Pref. | Kuro-sansyu-mikage | 4.91 2.89 24 2.65 
i | Ehime-Pref. Oshima-ishi 4.45 rth G spa ae is 
8 Kagawa-Pref. Yoshima-ishi 4.86 3.00 1G, Ped ith 
9 ” | Agi-ishi 4.56 2.78 21 2.71 
105 ” | Shokai 5.40 3.05 oer 2.70 
11 | Okayama-Pref. | Kitaki 4.83 2.96 20 2.68 
12 | ” Sabi-ishi 5.02 3.10 Aly 2.02 
13 y a 5.30 515005 . 26 2.94 
14 ” | Inushima-ishi 5.18 2.88 28 2.87 
15 Y Aoki-ishi 4.09 2.51 ay4{0) rea if! 
16 " Kitaki-ishi 5.71 3.55 .19 2.70 
17 ” Mansei-ishi 4.77 2.81 24 2.70 
18 | y Muneage-ishi 4.62 2.85 .19 2.77 

(1) Time mark. (8 microseconds.) 

(2) 300 ke/see. X-cut crystal. 

ti (3) 150 ke/see. X-cut erystal. 
(4) 50 ke/see. BaTiO3;. 
The average values of Porsson’s ratio as 


Fig. 6. 


Longitudinal 


i 
iM 


Waves in a 
specimen excited by different frequencies. 


sandstone 


determined by our experiments agree with 
Bircn’s value (1938). The averages of the 
ratios of the longitudinal and transversal wave 
velocities are as follows: 


ratio 
Granite | 1.69 
Marble | 1.78 


Finally, no tendency has been found that 
the velocity of elastic waves depends on the 
ultrasonic frequency within our frequency 
range (50 KC-300 KC) and also within the ac- 
curacy of our experiments. One of the ex- 
perimental results which indicates this con- 
clusion is shown in Fig. 6. This record was 
taken with a sandstone specimen. Although 
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Table II. Locality, Proper Name, Longitudinal and Transversal Wave Velocities, 
PorIsson’ 


s Ratio and Density of Marble. 
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ae ee | oy eee | ee ee hese) a Papen 
1 Yamaguehi-Pret. Akiyoshi Shiro | 5.62 a6. | 2.26 2.66 
| ” ” “Hakuehuzitu 5.97 3.32 28 2.73 
3 ” ” Midoriiriishi 5.34 2.96 .28 2.80 
4 ” ” Shirakumo 6.08 3.738 .20 hella} 
5 z ” Usugumo 5.48 2.92 .30 2.76 
6 ” ” Hakuyo 6.14 pe 31 2.76 
ul ” 7 | ” 5.09 2.80 .28 2.738 
8 " " Oka sarasa 5.89 2.96 33 alka) 
9 ” ” Kinmon 5.85 3.46 23 229 
10 " Y Kinteki 5.89 2.83 35 2014 
11 ” 7 ” 6.02 3.64 24 22 
12 ” ” Oka 3.76 2.02 ee 2.68 
1S ” ” ” (yeaalk Bt allies PAL 2.66 
14 ” Omine Hakuyo 5.25 3.13 24 2.83 
15 ” ” Shisuiseki 6.15 3.84 18 2.66 
16 Ua Vi Uzura 6.05 3.48 .26 2.712 
17 " Ota Sarasa 6.37 Bel “OL ZO) 
18 ” ” Okibuse 5.92 3.64 .19 2.80 
19 Tokushima-Pref. Tokushima Kamo 6.48 3.92 sya 2.76 
20 | ” ” | Kinkamo 6.34 3.43 29 2.76 
21 ” ” Aka-kamo 5.98 3.04 . 38 2.76 
22 ” Akaishi = Sato-shima 5.80 3.07 81 ee 
23 ” ” _ Kinkaseki 6.10 3.70 127 2370 
eit Kochi-Pret. Kochi - Ao-keiryu 5.41 3.15 2d 2.92 
25 ” ” | Keiryu 5.67 3.10 .29 2.82 
26 | Gifu-Pref. Akasaka Ki-sarasa 6.94 3.86 .28 2.94 
27 v ” Aka-sarasa 6.16 3.76 spalh De) 
28 ” ” rt 6.10 8.70 avail PARKS) 
29 " y | Kuro-sarasa 6.01 3.50 24. 2.62 
30 | ” ” Kuro 6.19 3.75 220 2.70 
31 | Ibaragi-Pref. Mayumiyama Shimanezumi 6.23 3.38 nay) AH 
BA | Okayama-Pref. Naruha _ Aiko | 5.12 3.18 .28 ane 
33 | Wakayama-Pref. Yuasa | Kamo | 6.02 3.10 382 2273 
34 | Yamanashi-Pref. | Shihanseki | 5.14 ig 80 rae 
35 Ryukyu Islands ” | 65.88 3.00 .28 2.89 
36 | Korea, Heizyo | Akagasumi | 6.00 3.30 .28 ao40 
aT) | Ttaly | | oan 3.51 26 2.66 
38 ” 5.90 3.28 28 2.58 
39 ” 5.70 8.20 PA 2.75 
5.53 299 was) 2.59 
oy ‘ ¢ 2.92 28 2.62 
Al a 5.25 2.92 ay 
42 | Va 4.68 2.438 5 a8 Fe 
43 Va Rok 2.75 330% ae 
le" 6.62 3.31 88 2.64 
ie : | 5.20 2.58 BB ao 
5.74 2.84 34. Paths 
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three different frequencies were used as 
indicated in the figure, no difference is seen 
in the time between the imput and arrival of 
waves. 

So far, the results for granite and marble 
have been described. We have made similar 
experiments for about 150 specimens of other 
kinds of rock. The results of these experi- 
ments will be reported in the near future. 
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Fluid Motions near the Core Boundary and the Irregular Variations 


of the core into mantle is studied. 


field near the boundary of the core are large enough to produce the observed variations 


in the rate of rotation of the mantle. 


the boundary of the core shows that the equations of the geostrophic wind in meteoro- 


logy can be applied to the layer. 


Bo 
in the Earth’s Rotation. 
By 
Hitoshi Takrucni 
Geophysical Institute, Faculty of Science, Tokyo University, Tokyo. 
and 
W.M. Exsasser 
Department of Physics, University of Utah, Salt Lake Giy,, Uta, Ur saae 
Abstract 
The transmission of mechanical forces of the magnetic field through the top layers 
It is shown that moderate fluctuations of the toroidal 
A simple theoretical analysis of the motions near 
have been developed to account for the 


§1. It has been suspected for some time 
that the observed irregular fluctuations in 
the earth’s rate of rotation (D. Brouwer, 
1952) can be explained by an exchange of 
angular momentum between the earth’s man- 
tle and core. W. Munk and R. Revett (1952, 
see also E. H. Vestine 1953) have studied all 
the available geophysical evidences and found 
that this explanation is by far the most 
plausible one. E.C. Butuarp and collabora- 
tors (1950) have shown that the observed 
westerly drift of the geomagnetic secular 
variation requires the existence of an electro- 
magnetic couple between the core and mantle. 
In Bullard’s model the earth is represented 
by three concentric solid spheres in mutual 
electrical contact and rotating relative to each 
other. In this paper, we go a step further 
by using the actual hydrodynamic equations 
for the region of the core near its boundary, 
and thereby investigate the transmission 
of electromagnetic torques from the core to 
the mantle. Simple numerical estimates show 
that the mechanical forces exerted by the 
magnetic field in the neighbourhood of the 
core boundary are at least as large as any 
other force acting in the fluid. Hence this 
problem involves hydromagnetic theory and 
is closely related to the dynamo models which 


earth’s magnetic field (W.M. EusassEr, 1950; 
H. Taxeucui and Y. SxHimazu, 1952, 53 and 
54; E. C. Butuarp, 1949 and 1954). 

As the following treatment shows, the 
motion near the boundary will break up into 
large-scale eddies, similar to the major per- 
turbations of the atmosphere. This is due, 
mainly, to the deflecting action of the Coriolis 
force, which leads to instabilities and eddy 
formation. In a previous paper (W. M. 
Extsasser and H. Taxeucut, 1955), we have 
approached this problem by using averages 
over circles of latitude and Jumping the 
effects of the eddies together in a mechani- 
cal and a magnetic eddy viscosity, which 
both are very large. It seems desirable to 
show how far we can go in a quantitative 
analysis without making this radical assump- 
tion. In the present paper, the frictional 
terms are considered to correspond to mole- 
cular viscosities and are small compared with 
other terms, especially the Coriolis and pres- 
sure terms in the equations of motion. This 
does not lead to a specification of the mean 
flow in the boundary layer, as in the preceding 
paper ; nevertheless it is possible to show that 
electro-magnetic torques transmitted from 
the core to the mantle are of the right order 
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of magnitude to account for the observed 
variations in the earth’s rate of rotation. 

§2. Let us take a local Cartesian coordi- 
nate system in Fig. 1, which has its origin 
at a point on the surface of the earth’s core. 
The equations to be solved are 


z: up 


y: north 


Bigs ak 
OH 1 y?H+ curl (Vx HA) (2.1) 
Ot 4rp0 
div H=0 (2.2) 
div V=0 (2.3) 
Be graa( ? +82 )—2(o x V) 
ot 0 / 
+ (curl Hx H)+yp?2V, (2.4) 
4np 


where H is the magnetic field, “4 the perme- 
ability, o the electrical conductivity, all in 
electro-magnetic units. V is the velocity 
field, w the rotation vector, p the pressure, 
p the density. Expecting that the boundary 
layer, where the westward drift of the earth’s 
magnetic field is taking place, is thin com- 
pared with the lateral dimensions concerned, 
we shall assume . 


(2.5) 


This will essentially be the only assumption 
required in this paper. Using this assumption 
in (2.2) and (2.3), we find that H, and V, are 
independent of z: 


(2.6) 


Since the velocity V, is relative to the earth’s 


mantle, which is in the region z>0 in Fig. 
1, and since V,=0 at z=0, we may take 

V.=constant=0 (C4500) 
in the boundary layer. Inserting (2.5)—(2.7) 
into (2.1), we get 


0H, 1 @H:z OV: 


alae H. (2.8) 
OL Eo Oz bi 0z 

Oy _ i Oak EF OVy (2.9) 
Ot. A4nyo 027 a Oz 

OH: _¢ (2.10) 
at 


By (2.6) and (2.10) is meant that Hz ised 
constant, independent of both (z,, z) and 7@. 
Using the preceding relations, we find from 
(2.4) 


OV: Op Wek Olaks 
+up?Ve, (2.11) 

OVy _ Op" _ nH. OH, 

gpentgy een lort east 
+up?Vy, (2.12) 


x CHOPS wiegt motpe rg c 
0=— + 2ey Ve ae (H,?+H,?) 


(2.13) 


where pp’=p+ gz is the pressure deviation 
from the corresponding hydrostatic pressure. 


The reason why we retain the terms — op 


Ox 
ap >. : ; 
ay in the above equations is that 
Bi 


although @ and 
Ox 


andar 


3 are assumed to be small, 
Ly 


p might be large and might keep —°P and 
v 
Op 


Oy. finite. 
§3. In the stationary state, we get from 
(2.8) and (2.9) 
OVe 1 = lho Only 
dz A4rnpoH, 02 ’ 
OVe ot Orde 
02 4npnoH. O27 eR 


Integrating (3.1) with respect to z, and putt- 
ing (Vz)z-0=(Vy)z-p=0, we find 
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Me el [off_( Atte) | 
"Anno H,| Oz OZ FG 


Ene 1 Pais 
4npoH,| Oz Oz Me 


In our present problem, H, may be considered 
the z component of the dipole field in the 
boundary layer. Since H.<0 (or >0) in the 
northern (or southern) hemisphere, the first 
equation in (3.2) will give the following con- 
dition for the stationary westward drift 
Vo<0: 


* 


(3.2) 


OH, (0H, 
Oz ( 0z de 
<_ 0: northern hemisphere 
> 0: southern hemisphere 
(3.3) 
The distribution of Hz as shown in Fig. 2 
will satisfy the above condition. Thus we 
have the distributions of the fields shown 
schematically in Fig. 3. These are nothing 


be a 
northern hemisphere southern hemisphere 


D 


Ee ————— 


RigMi2: 


but the distributions which have been sup-: 
posed hitherto to exist in the boundary layer 
of the earth’s core. The first equation in 
(3.2) will also give a quantitative result of 
the following type 


A,~4r p06 Ve DE, (3.4) 


where D is the thickness of the boundary 
layer and H, is the toroidal magnetic field 
at the bottom of the boundary layer (see 
Fig. 2). (3.4) is of the same type as obtained 
by Butuarp (E. C. Butuarp, 1949, especially 
his equations (1) and (2)). In Butuarp’s paper, 
H, and H, are written as H; and A, res- 
pectively. Putting w=1, o=3x10-§, V,= 
0.03 and H,=4 (all in c. g.s.e.m.u) in (3.4), 
we find 


H,~0.5D, (3.5) 


where D is meassured in Km. (3.5) will 
give H,=5, 25\and 5) gaviss tor D=I10) 50 
and 100 Km., respectively. These are plau- 
sible values for the boundary layer. Insert- 
ing the above results into (2.12), we see that 
the second term on the right-hand side 2a, 
V.~10-4x3x107?=3x10-® is about 10 times 
the third term pea 2S (H, was taken 
OPE — Xo) 
to be 4 gauss). Although we have no ex- 
perimental evidence on the magnitude of V,, 
we may assume that it is of the same order 
of magnitude as V,. Thus, in order that 


the fourth term aie on the right-hand side 


of (2.12) is of the same order of magnitude 
as the third term, we must have 
-»H2HyD 

Wes 


With D=10Km., V,=0.03 cm/sec and H.= Hy 
=4 gauss, (3.6) will give y~5x10°, which is 
much larger than hitherto supposed. For 
smaller values of v, the fourth term will 
become very small compared even with the 
third term. Thus, in the stationary state, it 
is very likely that the second term balances 
the first term and we have 


V~10 (3.6) 


re 2? 200,Va=0 (3.7) 


0 oy 
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Similary, we get from (2.11) and (2.13) 


m0 5 V0, 
p 0x 
Op” em 0 HH? H. OE 
_ +20, Va 8x0 aes Peat y’) a 
ee cal} (3.7) 
p 02 


(3.7) is nothing but the equation for the 
geostrophic wind in meteorology, and (3.7) 
and (3.2) are the fundamental equations in 
the stationary state. 
§4. On eliminating H,, Vz and V, from 
(2.8), (2.9), (2.11) and (2.12), we find 
0 ; 
2s L H,z=(p’) 


where (f’) is a functiou of p’ and 


z=| a : oO? rm oO? ae of FE 
7 ee az /\at ” dz 02" 


(4.1) 


0 0? \3 
02?|——Vm 5 
ee © > = sie 
ore cdupitely 
™ Arpo’ ie Arp 49) 


If we put w,=0 and v=ym=0 in (4.1) and 
(4.2), we have the following homogeneous 
equation for H, 


Oma OP 
—_—¢c?—— |\H,=0, , 
ee : Oz ) aes 
which will give 
fal. co ein(z— et) (4.5) 


This is nothing but the equation for a mag- 
netohydrodynamical wave first obtained by 
H. Aurven (1950). When w,, », ym2<0 the 
homogeneous equation of H, in (4.1) will 
give 


tale co et it ting (4.6) 
[—q+yn’)(—@+ mn) +en?? 
+407(—q+vmn’)?=0 (4.7) 


This is the expression for the characteristic 
time t=(real part of q) in which the stationary 
state will be reestablished after having been 
disturbed by a disturbance of wave length 
ya2™ 


Putting the values in the neigh- 


bourhood of (3.5), p=10 and 2=4D (rather 
tentatively) into (4.7), we find the values of 
the longest characteristic times t in Table 1. 


19) ENG H, : gauss 7 3 year 
10 5 0.5 
50 25 12 
100 50 50 
Table 1. 
§5. In this section, using the results ob- | 


tained above, we shall study the problem on © 
the irregularities in the rotation of the earth, 
which has recently been discussed in con- 
nection with the westward drift of the earth’s 
magnetic field (W. Munk and R. REVELLE, 
1952 and E.H. Vesting, 1953). Since the 
boundary layer is very thin, we have effec- 
tively two layers, the mantle and the core. 
Thus denoting the moments of inertia and 
angular velocities of the mantle and the core 
by In, I; Om and w,., respectively, we have 


Im4@m+1e4a0.=0 2 (5.1) 


where 4wm and 4w, are the deviations of the 
angular velocities from their stationary values. 
Since Im (== 7.2 x 10“) > J, (= 8.6x 10%), (5.1) 
will give |4w.|>>|4om|. Thus, we may 
write 


£0, = ere = 3 de A(M—m) 5 
m m 


(5.2) 


where 4(@¢—m) is the variation of angular 
velocity of the core relative to the mantle, 
which may be observed as the variation of 
the westward drift of the earth’s magnetic 
field. On the other hand, from the relation 
mT =2z, where T=86400 is the length of a 
day in sec., we find 


T? . Lae, 

AL= : i. Om == —-* == 
rmlaha = 5S edhe Wm) 

or 
eT: it Tie Lasees 
A(@e— Wm) fear we radian 
millisec. 

=7.9 degree , (5.3) 


decade 
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which agree quite well with the formula ob- 
tained by Munx and REvELLE using a little 
different model than ours. In their formula, 
msec. 
8.2 deg stands for our 7.9 in the above 
decade 

formula, the difference of which is only 
trivial. The stationary value of wo—am is 
—0.18°/year (E.C. Buttarp, 1950). Putting 
A(@e—@m) in (5.3) equal to 0.18°/decade, that 
is, one tenth of the stationary value, we 
have 4T=1.4 msec. This means that when 
the westward drift of the earth’s magnetic 
field is changed by one tenth of its stationary 
value, the length of a day will be changed 
by 1.4 millisecond. According to D. Brouwer 
(1950) (see E. H. Vestine, 1953, Fig. 5), the 
length of a day (a year) was increased 
linearly by 3.8 msec. (1.0 sec) from 1890 to 
1915 and decreased by almost the same 
amount from 1910 to 1930. These changes 
correspond to those in the westward drift 
of the magnetic field by 0.3 to 0.4 of its 
stationary value. These are plausible varia- 
tions to occurr during the above periods. 
With the symbol in (5.3) the above results 
may be written as 


AT=at, where @a==———_ (5.4) 
20 year 
or 
BAe 2 at 50x 10-2 (5.5) 
dt be 


Multiplying with the moment of inertia In= 
7.2x 1044, we get the torque acting on the 
mantle 

(6) 


fo oo 10 C2 .S. 


dt 
We shall now show that moderate fluctua- 
tions of the toroida! fleld near the boundary 
of the core are adequate to produce the 
electro-magnetic couple of the above order 
of magnitude. If the earth’s mantle were an 
insulator, the toroidal field there would vanish 
and there would be no electromagnetic couple 
acting on the mantle. It is, however, a kind 
of mathematical idealization to assume that 
the mantle is an insulator, For the con- 
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ducting mantle, the electrical conductivity of 
which may be much smaller than that of 
the core, a fraction of the toroidal field in 
the core will penetrate into the mantle in 
the way schematically shown in Fig. 4. In 
spherical coordinates (7, 0, ¢), the polar axis 
of which coincides with the rotation axis of 
the earth, the ¢ component of the electro- 
magnetic force Fy and the electro-magnetic 
couple 7’ acting on the mantle are given by 


--> 


core boundary 


Se 
es) 


mantle 


TAR 


[Yee eile 
Ho 0 , 0 He 
Fy= aD igd cere 3p Gye 
Haas See Gis ora a 
FH, OH Or, ; 
7} ——" _— (rH ae Ge’! 
Beles ar Se )| ce 


R+AR (2 (24 
r=| | \ Fe-Psintfdrdidb, (6.8) 
0 


R 0 


respectively. In (5.8), the integration is carried ~ 


out in the region between r=R (the core 
boundary) and rv=R+4R, where the toroidal 
field Hy is appreciable. In the above region, 


we have 
oe! 1 Pasay 
Or~ xr’ rsind’ rsinddb R 
Hey rewed eee (5.9) 
or 
ete wee Hon sin? 0d dé 
0 JO 


(5.10) 


In order to get a concrete result on the 
net torque on the mantle, we speciallize the 
poloidal and toroidal field to be a pure dipole 
and a pure quadrupole, respectively (W. M. 
Eusasser, 1946 and 47), say 


Rx 
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aP,(cos 6) 


i; = b, P,(cos 6), Ag a be F- 00 


(5.11) 


where P, and P, are Legendre functions. 
Inserting (5.11) into (5.10), we find. 


Tis : Rbybe (5.12) 
Equating (5.12) to (5.6) and putting R=3.5 
x10%cm, and b,=4 gauss, we have 

b,=90.05 gauss , (5.13) 


which is small enough to be a fraction of 
the toroidal field in the earth’s core. 
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A Simple Method for Calculating the Deflection of the Vertical 
From Gravity Anomalies, with its Applications to 16 Selected 
Stations in U. S. A. 


By 


Chuji Tsusor* and Akio Hayatu 


Department of Geophysics Faculty of Science, Tokyo University 
Tokyo, Japan 


§1. To calculate the —- and 7-components 
, of the deflection of the vertical at a point 
situated on the earth’s surface from gravity 
anomalies which are known on the same 
surface is one of the classical problems in 
physical geodesy since the day of C.G. Sroxxs. 
The well-known formula which was derived 
by Sroxes for this purpose is theoretically 
beautiful but unfortunately it cannot be used 
in its original form at present because of our 
rather poor knowledge concerning the gravity 
anomaly distribution all over the world, which 
is essential for the application of that method. 
Among the many efforts which have since 
been made to make this method of Stoxrs 
competent for practical applications, the deri- 
vation of a useful formula by Vening MrInesz 
(1928) is to be particularly mentioned. In 
this formula, the effect of the gravity anomaly 
in an elementary compartment of the earth’s 
' ‘surface upon the deflection of the vertical at 
a certain point on the same surface is given 
in terms of the angular distance ¢ of the 
compartment from the point and of its azi- 
_muth @ as seen from the same point. The 
total deflection is obtained if such an ele- 
mentary effect is integrated over the entire 
surface of the earth, that is, from 0 to z with 
respect to # and from 0 to 2z with respect 
to @. 

Souuin’s table (1947) which was published 
later is an important contribution to this 
field of study. In this table are given numeri- 
cal values of the quantities needed if one 
wishes to apply Metnesz’s formula to actual 


problems. 


Making use of this table, Donald A. Ricr 
of U.S. Coast and Geodetic Survey (1952) 
calculated the &- and y-components of the de- 
flections of the vertical at 16 stations selected 
in U.S.A. using free-air anomalies known 
around them. Rick gave curves showing how 
the calculated & and v tend to their respective 
definite values according as the limit of inte- 
gration with respect to ¢ is successively in- 
creased. Looking at these interesting curves 
of Rick, it was noticed that the greatest dis- 
tance taken by Rick in the direction of ¢ was 
a little over 600km at most, which is onlya 
small fraction of the earth’s radius. So far 
as Rick’s calculations are concerned, there- 
fore, no appreciable part seems to have been 
played by the assumption of the sphericity 
of the earth which underlies Mrtnesz’s for- 
mula. It naturally occurred to the present 
writers that it will be interesting to see what 
will happen if the sphericity of the earth is 
ignored altogether in this problem. There 
can be two ways of approach to investigate 
this, say (A) and (B): the one (A) to solve 
the plane problem from the outset, while the 
other (B) to find the approximate expression 
of Metnesz’s formula in case when ¢ is made 
small. It might well be expected that the 
two will turn to the same formula after all 
and that this formula when applied to the 
actual problems in U.S. A. will give results 
that differ little from those obtained by RIcE 
with much higher approximations. 

§2. First, we shall follow the way (A). 
Taking the cylindrical co-ordinate system 
with its origin at the point at which the 


* The preliminary report was published in Proc. Japan Acad., 30 (1964), 461. 
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deflection is to be calculated, the gravity 
potential U (anomalous part only) can be ex- 
pressed by a Fourt=r-BEssEL series such as 


U=>) i ao In(kr)exp (—kz)dk-cosnd 
mw JO 


in =|" io Inher) exp(—ka)dk-sinnd ....(1) 
n Jo 


using customary notations. From (1), the 
gravity anomaly 4g at any 7, 9, and gz is 
given “ 


ag=— a Ds i Be, Jn(kNexp(—kz)dk+cosnd 


+E\" Bi, Jn(enexp(—ka)dh-sinn® . .. :.(2) 
n Jo 


Since z=0 at the earth’s surface, (2) reduces 
into 


4g= x \° B¢, Jn(kr)dk- cos nd 
= i Bi Julkdk-sin nd 
= ie Bé, Joker dk 
+ |" Bian(kridt-coso 
* \ BS, Jo(kr)\dk- cos 20 +--+ - 
- \, Bs J(krdk-sind 


aM Bs JARr)dk-sin26-+.-.. 
0 


=C)(r)+Ci(n)cos 0+C.(r)cos20+---- 
+S,(7)sin@+S,.(r)sin20+.-... 


where Cirn=\" 25 Io(krdk , 
0 

cin=|" Be. J,(kr)dk , 
0 


sins |" Be jkndk , 
1 


On the other hand, the & and y-components 
of the deflection of the vertical at the origin 
are given if y in the following expressions, 
(4) and (5), is made to tend to 0. 


= (or be 


Chuji TsuBo1 and Akio HAYATU 


§ 
=U" Bi.nien dtr |"Bian Erde, 
0 
eed 
ah or )e 
= Or 6=2/2 . 
1 | 
=7\\2 Beh (edhe) ey | 
& 
If, in doing this, the following mathematical 
relations concerning BrsseL functions are 
noticed, 
; 1 / / 
jh Ongr IO=) (0)=----=0, .- (6) t. 
the following expressions are obtained, 
Lice Ss 
=— bye 
E rap Be, dk, 7) 
pee eh nC) 
2g 4 1k 7 
By, and Bs, being the constants which already 
appeared in the expression (3). The next 


task then is to evaluate the integrals in (7) 
and (8). If the following mathematical rela- 
tion again concerning a Brsset function is 


noticed, 
ne dy=T; 
0 ug 


the integration with respect to k can be re- 
placed by that with respect to 7, thus 


\Bidb— \ \. Be, WD) ae ay 
0 0 


« 9) 


-|" GW) ay, (10) 
0 Ve 
and similarly 
Byak=| Slope xii: 
0 Cay 
The final expressions for — and y are there- 
fore simply 
Ci(r) 
E= aa : dr, pepe > 
Sir) 
=26\, : dr beh) 


These expressions state that what are needed 
for our purpose are only C,(7) and Si(7) which 
are the first order coefficients when the 
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gravity anomalies along a circle r=r drawn 
around the origin are expressed by the 
Fourier series of azimuth. So our method 
for calculating & and y practically consists 
af - 

1) drawing a series of concentric circles 
around the point for which € and 7 are to 
be calculated, 

2) reading the values of 4g along each of 
the circles at a reasonably many number of 
points, 

3) calculating the first order coefficients 
C\(7) and S,(v) for expressing 4g along each 
of the circles by the Fourier series of azi- 
muth, 

4) dividing C,(7) and S,(r) by correspond- 
ing 7, 

5) integrating C,(r)/r and S,(7)/r with res- 
pect to 7 from 0 to oo and lastly, 

6) dividing the integrals by wen 
Ci(7r) e “a 27) ay may 


first appear to diverge ate 7 = 0, ae r is 
contained in the denominators of the inte- 
grands. But this is not actually the case, 
because C,(r) and S,(7) themselves tend to 0 
also as the first power of 7 if ry tends to 0. 
If v is small, C,(7) is a linear function of 7 
and may be written as 


The integrals pe 


C\(r)=ar, VS 7 Te cla) 
where F is the limit up to which C,(7) may 
- be regarded to change linearly as ar. The 


integral from 0 to co may be split into two 
parts: the one from 0 to R and the other 


from R to o, and may be written as 


C Cur) 
Fe eae Ti 


Ls ate C7) sf nt Ci(7) dr| 
20 fe r = Of 
ae \" Co) ar] 
28 | nT 
1 f °- Ci(7) 
S75 R SS eatiliar 
ag + |" - 
The value of & as given by (15) is evidently 
finite. The splitting of the integrals into the 


two parts, however, need not be necessarily 
exactly at the very limit up to which C,(7) 


(15) 


AT 


changes linearly with rv. If, by any method, 
E(R’) is known to be the total contribution 
of dg within r=R’, — may be written as 


E=E(R) +, al med r, HG, 
where R’ may take any value. Also 
2= WR) +5— pele tee Bul) ap Pal 
R’ 


§3. Now in order to investigate its feasi- 
bility, our new method will be applied to 
calculate the values of € and 7 at the same 
16 points in U.S. A. for which Rickr’s values 
have been known. The comparison of the 
values by Rick and by the present writers 
will make it clear how far the new method 
can safely be used. 

On the map of 4g published by Ricz, concen- 
tric circles were drawn at every 20km distance 
interval around the point at which € and 7 
are to be calculated. C,(7) and S,(7) were 
calculated from 24 equidistant readings of 4g) 
along each of the circles around the origin 
by the ordinary harmonic analysis method. 
As to the values of &(20) and 7(20), however, no 
accurate estimations have been possible for us 
for the map published by Rick was_ un- 
fortunately of too large a scale to admit them. 
vee could be done was only to calculate 

AN Cir) dr and Ale ree 
250 7 2g) 7 
(16) and (17), and to the calculated values to 
add the values of & (20.09) and 7 (20.09) ob- 
tained by Rick. At the writers’ request, the 
values of £ (20.09) and 7 (20.09) were kindly 
supplied by Rice. He sent also the values 
for larger distances for reference. 

In order to demonstrate the procedures 
followed in our calculations, the results for 
Meades Ranch, one of the 16 stations, are 
given in Table I as an example. 

Since the integrands in (16) and (17) contain 
y in the demominators and C,(7) and S,(7) do 
not increase systematically at large 7, the 
integrals decrease with vr. Owing to this 


GM ay an a\ Si) dy 
20 ¥r 


do not change much at e. y, so that the 


dr, according to 


reason, the values of (i ‘ 
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integrals from 0 to co can be approximated 
by those from 0 to a reasonably large 7. 
If &(7) and 7(r) have the following meaning : 


E(r)=£(20)-+ » \, OO) yp 
jo 7 

hr) =4(20) + mar Sit) gen a pd) 
Zjo 7 


Table II will show the values of &(7) and 7(7) 
found by our calculations for all the 16 points. 
The values by Rick are also entered in the 
table for the sake of comparison and are 
printed in gothic type. 

TABLE I 

Meades Ranch 

Component 


| Yd (& | 
(km) Cy(7) Ci(r)/r \ Cur) dr E 
| (magl) | 10-9 | 10-a (second) (second) 
= Si - a | mere Fe al. pai 
0 be vel 
207 0.33 1 "0.17 — —| -0.28* 
40 6.42 1.61 Abr 0.19 | -0.09 
60 6.40 1.07 4.44 0.47 0.19 
80 2.26 0.28 5.79 0.61 | 0.33 
100] 3.81 | 0.38 6.46.|. 0:68 | 0.40 
120] 112] 0.09| 6.93 |~ 0.78} ~0.45 
140 | -2.86 | -0.20 | 6.82 0.72 | 0.44 
: | 
160 | -~2.76 —0.17 6.44. 0.68 0.40 
180 | -4.08 -0.23 6.05 0.64 | 0.36 
200 | ~7.34 -0.37 | 5.45 0.574)\an 0:29 
220 | -10.78 —0.49 | 4.60 0.48 | 0.20 


* RICE’s value 


n-Component 


OS ee eee ee 


* | Sir) | Sor | Par) 4 
em fe 

0 v 
20| 5.48] -2.74 a = be 0) D6R 
40 | -8.66| 0.92 865g 369 34 
60 | 0.11 0.02 4.55 | 0.48 | -1.44 
80| 4.29] 0.54 =4.00 | 0,421 21,98 
100.| 7.49| 0.75| -2:71| ~0198 | 1294 
120! 0.50! 0.04 “1,92 \npe20 11S 16 
140 | 6.16] -0.44 ~2.82 | -0.24] -1.20 

160 -4.11 | -0.26 eB OH ee SO Veet 
I> 4. ; 02 | 0. ~1.28 
180 | -4.34] 0.24 8.51 | 0.7 heet oe 
200, -4.39 | -0.22 -8.97 | ~0.42 | -1.38 
220. 1.37 | ~0.06 2A, 26. 0204501. 


* RICEH’s value 


The values of &(7) and 7(7) are plotted in 


Fig. 1 against 7 up to which integrations 
were extended, It is clearly seen that the 
values of Rice and ours agree satisfactorily 
well, At a few points, Little Rock Longitude 
for instance, however, there is some systema- 
tic difference between the two which amounts 
to 0.3”. The origin of such a difference is 
not yet clear. 


§4. Now the second way of approach (B) 
will be taken to see what will happen if ¢ 
is made small in Mertnesz’s formula. The 
formulas as given by Metnesz are as follows: 


fa ei df (¢) .: 
pas yt Ag (v0) =~ sin¢cosOd¢dé, 
ang \ \ isis dy 
oye (20) 
1= 5), | 4g(¢0) GAS) srgcinbde do, 
27g 0 Jo d¢ 
ae 4 
7 i 1 Meades Ranch 
oh ee ee tate 
e 
DNerenyeens be e 
e 
e 
e 
e 
pa) 
27 4 Twin 
| Pe e ® 
AY ; MO 300 400 500 600 mitts 
e 
a ar te 
e 
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of ah = 608 Hid Joke cos 25 sin ¢ 
dp 2 2sin 29/2 a) 


+3 sin ¢ log (sin § +sin 
2 2 
-* noe _ cot 4 cos | wpe heed) 


in which ¢ is the angular distance between 
a compartment on the earth’s sruface where 
4g is known and the point at which the 
deflection of the vertical is to be calculated. 
@ is the azimuth of the direction that connect 
the two. 

As already stated by Sonuins (1947), if ¢ 
is small, (22) can be written approximately 
as 


df il Sy 
a —— : ae cul VX 
dv y? 20 2) 
Therefore 
6+ (2402 1 
E=5, | | Ag(o0)*cosddgd0, ....(24) 
2xg)Jo Jo g 
and — 
1 27 (1 iL s 
15,4 4g(¢0) —sin@d¢dé. 5 03(25) 
ang 0 Jo ¢ 
But since 
22 
=| 4g(¢0)cos 0 d0=C\(¢), (26) 
0 
and 
1 ee : 
= Ag(o0)sin 0 d0=S\0), (27) 
0 
& and 7 are given by 
pal! Cx (Y) 
E= AK 5 dg, x tee) 
and 
Si(¢) 
=3,\. F dd, eo ees) 


If, by the reason already stated, the integra- 
tion with respect to % is extended up to ¢ 
instead of up to z, (28) and (29) reduce to 


_1(* G@) 
g =s,\, 3 oe, (30) 

a . Sig) 
2g \ : dt. ..(31) 
On the other hand, if the earth’s radius is a, 


av is our vr. Thus (30) and (31) may be 


written as 
Ee =k Sle rs (32) 
a a 
val Sur) dr, . (33) 
~ Og r * 


which are nothing but the expressions already 
obtained in (12) and (13). 

Proceeding to the second approximation, 
if the term = in the expression (23) is 


retained, € and 7 are given 


1 
es 2, ly 3A 
=5 Ae, +57 )ar, 168 
y 
7 50\, SiO( 5 vt) 
=32\,5 7) ( +5,)4r- (35) 


The following is the summary of what have 
been sta‘ed, 
a) in the first approximation, 
f= I hes dr, 
22g x 
Si(7) 
dr, 
=|, F 
b) in the second approximation, 


1 eran db a 
=—. GE — —— 
E = (7) ( - +5, )ér > 


1 {= ie 
eos \ (45 jar, 


c) in the third approximation. 
df) 
ss ¢ 
p=) \ aw 
1 
ale Syl pees a sin ¢ dy 
& 
The three curves in a 2 show the values 
of 


i (ea 
—| — do 
=) eile ¢ $ 


1ty/1 3 
ete € 4 
save aa 5 dt, 


= ht) 


veel} 


sin¢d?¢, 


y= 


as plotted against ¢ or r. These are &s 


Ware 
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¢ 


c) J siny E ay 


rin km 
200 400 600 800 1000 
Oo rae See aS ee St $$ 
12 woes oy 6. 7. BB 9) “410 
Y in degree 
Hig: 2. 


when C, or S,; are dropped out in (36), (37) 
and (38). The three curves are very much 
similar with each other and it will be under- 
stood that even when C,(7) and S,(7) do not 
depend on 7, the error in the values of & and 
” introduced by our first approximation is 
very small and may safely be ignored. 


Thus the conclusion is that the formulas 
(12) and (13) proposed by the present writers 
are very satisfactory ones and are widely 
recommended owing to simplicity of use and 
reasonable accuracy which can be attained 
by them. But if one wishes to get more 
accurate results, the formulas (34) and (35) 
are to be recommended. These formulas give 
almost the same values as the exact ones do, 
provided, of course, C\(v) and S,(7) do not 
increase much as 7 increases. 

In conclusion, the present writers wish to 
thank Dr. Donald A. Rice of U S. Coast and 
Geodetic Survey for kindly placing his 
unpublished values of & and vy at the 16 
American stations at the writers’ disposal. 
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Equation of State for the Earth’s Mantle Based Upon a Theory 
of Pressure- and Temperature-Dependences of 
Elasticity of Solid 


By 


Yasuo SHIMAZU 


Institute of Earth’s Science, Nagoya University, Nagoya. 


The pressure at the center of the earth is 
generally accepted to be of the order of three 
millon atmospheres. The temperature there 


is believed to be several thousand degrees. 


In order to understand the physical state of 
materials at such a high pressure and tem- 
perature, experimental studies in high pres- 
sure and temperature physics are highly 
useful. But the highest pressure that could 
be produced so far in the laboratory is about 
1.5x10° atmospheres, which corresponds to 
the pressure at 300km depth of the earth, 
only 5 per cent of its radius. The highest 
temperature that could be produced experi- 
mentally under high pressures is only a 
thousand degrees. It must be admitted there- 
fore that these man-made pressure and temp- 
erature conditions are as yet rather far from 
those prevailing in the central parts of the 
earth. Thus theoretical considerations are 
naturally needed particulary concerning with 
the equation of state of materials of which 
the earth is composed. 

The materials of which the earth is com- 
posed are in a state of finite initial defor- 
mations caused by the high pressure to which 
From the propagation 
of seismic waves within the earth and the 
periodic deformation of the earth by tidal 
forces, we know that the mantle of the earth 
is an elastic solid. Seismic waves are pro- 
pagated through this elastic solid which is 
under finite initial deformations. Since the 
velocities of seismic waves are known at any 
depth of the earth from observations, we can 
estimate from them the order of initial defor- 
mations, if the behavior of elastic waves pro- 


_ pagated through the finitely deformed solid 


be known. 

Let us first try to extend the theory of finite 
elastic strains due to F.D. MurnacHan, L. 
Brinvtouy, R. Karrus, R. Fijrvs, etc. to obtain 
exact stress-strain relations *(1, IJ). In our 
theory, the strain under considerations 
must be perfectly reversible. We shall 
study the effect of initial dilatation and shear 
upon the elastic constants of both isotropic 
sohd (1° Il) and cubic crystal (1 V)=.the 
case when hydrostatic compression is work- 
ing is also included in our calculations. 
Elastic constants are expressed in terms of 
the initial elastic constants and of the applied 
initial strains or stresses. For example, the 
bulk modulus K and rigidity G have the 
forms as follows: 


Isotropic solid 


eo ( (24e80= 1816 6 n) 
22a 


3 Ko 
Pe 
G=G)+=—— (3449 u+15m+0.57)+...., 
3Ko 
Ay, 
Cubic crystal 


K=Kot+ pg (4Cn t8Cn—2Dis—4 Dns 


OG 
-; Diss) are wey 
1b sea 
G= Grrtgge (Cut 2Cnt7 Cu Dyan) 


ae Bee) 


where Ky, G) are their respective initial values 
and P is the applied hydrostatic pressure. 
(l,m, m) and (Dy, ..--) are the elastic con- 
stants of the second order. 


SoG bs II) denotes Chapter II of Paper 1. The present paper is the summary of that paper. 
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Now We shall discuss the propagation of elas- 
tic waves through anisotropic solid under initial 
strains (1, 1V). Since anisotropy is introduced 
when initial strains are applied, velocities of 
the waves become dependent upon the direc- 
tion of propagation. Isotropy is kept only 
for hydrostatic strain. The initial (unstrained) 
coordinates will be denoted by (a, b,c) and 


the final (strained) coordinates (of the same 
particle) by (x, y, z). When uni-axial tension is 
applied along a-axis, the displacement is ex- 
pressed by (x—a=aa, y—b=— Bb, z—c=—B¢). 
In this case, the velocities of elastic waves along 
several special directions of propagation are 
as follows: 


TABLE I 
Siegen cos. || displacement of | 
of the wave | a particle corres- i | 
propagation ponding to the | Veloriy _ Remarks 
wave | x 
| —a—8y —0.67 
4=1 | 18 eae Vel dr Aa eee) | purely trans- 
| ie versal wave 
—2va-8p-—m-0.5n 
| +e! ) } ! 
| : | 
—8aA-61+2m 
m=0 Leas ls aes | Wa= Vol $2 ropes ony | quasi-longi- 
ls qs ay tudinal wave 
| Spe (= A—14p+12l+2m | quasi-trans- 
Wra+2p) )t | versal wave 
. - —a—-0.5m— 
FeV We= Vell +e( ak BL 
ls Y2 24 
jf ~2a—Su-—m 
ba a ( 2m )} 
== 7 rae r -3 ti 5 | 
ie Wea Wa. e(— nor in polarized 
ae wave 
2 5 | 
e w—0.5n | 
Se 2 
( 2h ) 1 
a0 reer ss. Eee iene | 7 
| Ww eo Ww. , me "i —4a-14p-—2m | . 1 ii 
l.=0 * ale ne . 2aA+2n) ) ' polarized 
wave 
+e(Sti DEEN | 
2A+2u ) 
Wo=Wo | 
q3=0 Similar to the case q.=0 by putting g, for q3- 


eee Eee eee 


In Table I, Vay Aten and Vay are 


the initial values of Pelocities, € and = are 
defined by (1+ a@)?=(1—2e)-}, d—pyr=(1-— 
2e’)-'. Only the first order terms of ¢ and 
é’ are shown in Table I, whereas our cal- 
culation includes the second orders. Putting 


é=ée’=—f, the results obtained above be- 


0W?= 


come in the case of hydrostatic compression: 
oW y=A+2u+f (164+20 w—181—4 m) 
+?(62.52+ 65 »—2071—55 m 


—3n), 


ut (32449 w+1.5m+0.5n) 
+f(15 4+27.5 w—2714-1.5 m 


+2.5n), 


r 
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0=po(1+2f)3! si weet) 
For the finite plane shear which can be 
denoted by (e—a=tb, y—b=r’a, z—c=0), the 


results similar to those in Table I are obtain. 
ed as follows: 


TABLE II 


: 5 | 
direction cos. of | displacement of 
the | Wave propa- | particle peceeupGReihe 
gation | to the wave co 
ae s | = == = = = a ee 
l3=0, Wie Vp 
ee {aaah Pet a 
| %3=1 
pores Win= Vs 
Qo=1 similar to the case qg,=1 by putting - for -' 
i3=1 = Vip ; ry . 
B= 1 ls =0 
Wit at 52) 
: =+1 | Wr Vs| 1( pitts | 
| =U XS (i mm 
| Wi= Vp[1+(= ee ES: = et 2) 002 | 
q3=0 | ly, lg: complicated Wu=Vs 
38u+0.5 
I3=1 Win= Vel 145 (— a Nr + 2") 0 | 


It can be seen that the effects of initia] 


strains upon velocities of waves are small. 


Variation in velocities will be 1 per cent at 
most when a linear extension of 0.5 per cent 
or a shear of 40’ is applied. The range of 
perfect elasticity of materials of which the 
upper part of the earth is composed is much 
less than 0.5 per cent for extension. If any 
variation in velocities be detected before or 
after great earthquakes as some-one insist 
the cause for it is not likely due to the 
change in working pressure, and therefore 
it is beyond of the present study. 

AS we see. in? (1)) (2), sor Tables “I, UII, 
second order elastic constants appear in our 
calculations. All high pressure experiments 
to date to determine the orders of magnitude 
of these constants have not been sufficiently 
accurate except for bulk modulus. The only 


way left to predict the order of magnitude 
of second order constants is to develop a 
theory of finite strain for crystats (1, V). 
For crystals with a given lattice structure 
and with a given force law between consti- 
tuent atoms, the potential energy isa known 
function of deformation (=strain invariants) 
and the elastic constants of any order can be 
expressed as the coefficients of each term in 
this function. The potential energy between 
two particles of the lattice will be assumed 
to be of the form 


9 (nau | ne) + a) | .. (4) 
n—m m\ r nm\ 1 
(2 >m) 
By thermodynamical considerations, we get 
Seals ..5) 


(a) _m-+? 
OP / p-o 3 
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a) can also be obtained from the 
OP/p-0 


theory of finite strains as 


) 42 3DutDist6Din 
oP 3 Cy t+2Cy. 


It is confirmed by BripGman’s experiments that 


a) for both single and polycrystals is 
0P/p-o 


about 4 for 24 elements and 4-5 for 
40 compounds. These results lead us to the 
conclusion 


6<m+n<9 ye 0 


from (5). Then we get 


Z 
2Diyit+ 4D» sae Dy»3 


4C,,+8C,, 
in (1). From the above considerations, we 
see that the second order term in bulk modulus 
is nearly equal to zero. Unlike for the bulk 
modulus, the second order term for rigidity 
cannot be obtained by such thermodynamical 
considerations. Calculating the lattice sums 
for face-centered cubic crystals, we get the 
second order term for rigidity as 


Ae Bee a) 


2M Dog 
0.40>6.- +2C,47C, de LIS = av 39) 
corresponding to 6<m+n<9. Of course, 


OK 
oP 
lattice sums, and this gives the same result 
as (5). Simple and body-centered cubic cry- 
stals are found to be mechanically unstable 
for such comparatively small values of 22 and 
n. Wesee that the contribution of the second 
order term to rigidity is much larger than 
that to bulk modulus. Although experimental 


can also be obtained by calculating the 


0 
values for ( a) are scattered, we may say 
P=0 


4<m+n<l2. Then we obtain 


Corey’ gy’ n 
1.5m-+0.5n 
0,40>Lbm+0.5n. 9 1 
MME Treicc | ae 


for isotropic solid by comparing (1) and (2). 


Since we get the equation of state of solid 
as 


7/3 6/37] 
P=5 (e) aft [on#2e: 
We 2/3 
+) (Dir + Dust 6Dis){1—(5) Lead 


mean 


from the theory of finite strains, we may 


say that 
Vo ee) 
ele) GF) | 


3Dint+Di3t6Dyie =A 
2(Ci,+2Ci2) 


pele) 


and m=2, n=4, when 


from (8). 

Equation of state for a simplified model of 
the earth’s mantle in accordance with the 
theory of finite strains will be discussed as 
follows (1, VI). If the hydrostatic equilibrium 


is satisfied within the earth, density-pressure - 


(p—P) relation is found from 
Ks 


oP ee ee 
el Fad W3=o= ae 
where K; denotes the adiabatic bulk modulus. 
The relation between the adiabatic modulus 
and the isothermal modulus in (1) can be 
obtained by thermodynamical considerations. 
(13) is the fundamental equation in BuLLEN’s 
density distribution and we adopt it here too. 
In this model, chemical homogeneity is also 
assumed. The physical state of a material 
can be described by four variables : pressure 
P, density p, temperature T, and chemical 
concentration’ ” (SHimazuU, Y: 1952). The 
equation of state is usually expressed by 
p=o(P,T,n). Although there are many 
studies on the constitution of the earth based 
upon the density distribution, we must note 
that observable quantities within the earth 
are Wr and Ws, but not 9. Weshall discuss 
directly the earth’s model based upon Ws 


....(13) 


and Ws. From (13) we get 
OKs 1 d¢ 
a1 ne 
oP Oo ine i) 


Using seismic data, we see that a is nearly 


. re tegOR i 
— ~_ n —-— 


ea 
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constant (=3.7) for the layer between 1,000— 
2,900km depths. The layer shallower than 
1,000km depth is not in this consideration 
on account of an apparent scattering in the 


OKs 


distribution of - AP Then we get m+n2~5 


for the deeper layer of the mantle. Consi- 
dering the effect of temperature gradient, 
m+n may become slightly larger, say, m+ 
n~6. Thus the effect of second order elastic 
constants of the constituent material may be 
sufficiently small in virtue of (8), and the 
equation of state can be expressed by (12). 
Assumption of chemical homogeneity is also 
valid for this layer. 


Even if the second order terms in ¢ are 
nearly epual to zero, those in Wr and Ws 
may be large (1, VI). Calculating initial 
strain (f)-depth relation and using observed 
Wr, Ws distribution, we can derive Wr, Ws 
—f relation. We see that Wp and Wy are 
linear functions of f. Comparing these re- 
sults with (3), we get 


A=1.01 2, 
18/+4m=9.46 1 


1.5 m-+0.5n=—4.76p. -'...(15) 


We see that the effect of second order term 
upon the shear wave is greater than upon 


the longitudinal wave. 


In the preceding discussion, we see that the 
pressure is the predominant variable in the 
deeper layer of the mantle. It has been 
proved that the shear wave cannot be pro- 


~ pagated through the core of the earth. This 


may lead to the fact that fusion takes place 
within the core. Melting of materials there 
is related to the temperature. We _ shall 
study, for the next step, the effect of tem- 
perature upon the elastic constants of solid. 
Our main object is to calculate the tempe- 
rature dependence of rigidity of a solid. 


Now we shall study the behavior of an iso- 
tropic solid using Desyr’s approximation (1, IX). 
We may safely say that the mantle of the earth 
is in a high temperature region, that is, the 
temperature there is far above the Drsyr’s 


characteristic temperature @>( =") of its 


constituent materials. v) in the above expres- 

sion is the Drsyn’s characteristic frequency. 

Then the free energy density (=the elastic 

potential energy density) is given by 
oe ii On 


p= hot V 7 


.. (16) 


where ¢) and @p are expanded by e(=—/ in(3)) 
using (3). GriineIssEn’s parameter y¢is given by 


d log vo TaR 
Ge saa: t e= 2 
g Piet We get e KV, because 
oe must vanish because of the assumption 


that solid is in the state of zero stress at 
any temperature under no external force. 
The mean value of the coefficient of thermal 
expansion is then found as a= 2 _ 3Rra 
IP GES 
which gives the Grtneissen’s equation of 


OK : ; 
tate. —— t 
state eo 12 is obtained from 


(or). =(or).*Uoe Joon), 


Finally we get 


0K R 
14333 = |= 1.9169 ee 
re=1.833, (57) A Gas 
On, Anat LO We nouce thatenasicnal 


function of volume while it is independent 
of temperature. Bulk modulus decreases 
linearly under rising temperature at high 
temperature region. 

Suppose an isotropic solid is placed in a 
rigid tube. Raising the temperature, the 
solid will be expanded in one direction. There- 
fore, thermal strain is of the unilateral dila- 
tation type. Calculations of velocities of 
elastic waves propagated through the solid 
with finite unilateral dilatation are given in 
Table I. Using these result, we obtain the temp- 
erature dependence of A+2y. As it is 
very tedious to calculate the temperature 
dependence of rigidity using the finite strain 


theory for shear, we obtain toa by com- 
OT /p 


bining the results for A+2 and K=4+ 3 tte 
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The result of the calculation is (1, IX) 


ee 2418 R ‘Gea = 203.07" 
2 P 


OF Vo’ \OT Vo 
ao + AS) 
for 4=y at T=0, and finally 

Ok R OKs R 
—— ) =—19.69 — = =—9.61—, 
Caray, Gr)2 oy, 
0G R 

—— ) = —23.27 —. ee uke 
(ar), Ca ne 


Our result shows that the rigidity decreases 
linearly and reduces to about one-half at the 
so-called melting point. There are many ex- 
perimental works on the temperature depen- 
dence of elasttc constants. Attention has 
specially been given to the decrease in rigidity 
under rising temperature in relation to the 
study of the fusion mechanism. It has been 
shown that rigidity reduces to about one- 
half at several degrees below the melting 
point beyond which it decreases rapidly to 
nearly zero. That rapid decrease near the 
melting point and its relation to the mecha- 
nism of fusion are beyond the scope of our 
study. 
From the equation 


aa)-—" Laple) k*Laxte) bar’ 
BA .(20) 


we find (1, X) that the effect of temperature 
upon K or G, for 1,000—2,900km depths, is 
only 5 per cent or 13 per cent respectively 
of the effect of pressure. In that calculation 


aT is assumed to be adiabatic. Adiabatic 


dr 
temperature gradient is obtained as shown 
later. If the melting point gradient obtained 
by R. J. Urren, which gives thc maximum 
temperature gradient, is inserted into ar , 
7 
the effect for K and G is nearly 16 per cent 
and 40 per cent respectively. However, when 
we examine the thermal state of the earth’s 
interior, it is likely that the thermal convec- 
tion exists within the mantle, and that the 
actual temperature gradient would be nearly 
equal to adiabatic. Thus we see that our 


model for 1,000—2,900km depths is. valid. 
It is not confirmed, however, that our model 
is valid for shallower mantle. Finally we 
find that the temperature at which the rigi. 
dity of a rock reduces to one-half is 6,000 
°K at the pressure of the core boundary, 
whereas it is 1,790°K at the surface. At 
the core boundary, the adiabatic gradient 
gives the temperature of 1,600—3,200°K cor- 
responding to the assumed temperature of 
1,000—2,000°K near the surface. 

If the mean atomic weight of the consti- 
tuent material is assumed, Drpyn’s charac- , 
teristic frequency ») and GRUNEISEN’S para- 
meter 7, at any depth can be estimated from 
seismic data. We take A=20 from the mean 
value of rocks. Further, the distributions 
of the coefficient of thermal expansion a and 


the adiabatic temperature gradient (from 


d log cen a 
dr 


C, £) can also be _ estimated. 
p , 

There is a comparatively rapid decrease of 
a with depth. Using the distribution of v9, 
we can estimate the distributions of some other 
physical quantities. Theremay be an increase 


in electrical conductivit «x ex a “ 
y ( p RT 


W=W (vo)) and a decrease in viscosity (« exp 


)) with the depth. Thermal conductivity 
(ccol/3y9) may increase moderately to reach 
2.3 times the surface value at the bottom of 
the mantle. 

Is it likely that the thermal convection 
exists within the mantle? To find one 
answer, we study energy transfer within the 
mantle (1, X). According to a theory of the 
origin of the geomagnetic field, (TakEucat, 
H. and Suimazu, Y: 1952-54), the heat 
flow from the core to the mantle is of the 
order of 2x10-® cal/cm? sec at the bottom 
of the mantle. The heat which is transported 
through the mantle by thermal conduction 
may be, at most, of the order of 10-® cal/ 
cm’ sec. Once the convection takes place, 
the amount of heat transfer may reach 3x 
10-° cal/em? sec. This result suggests that 
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the convection is sufficient™ but the conduc- 
tion is insufficient to transfer the heat energy 
generated in the core to the surface of the 
earth. Worthly to be noted is that the mean 
heat flow through the surface of the earth 
outward is of the order of 1x10-° cal/cm? 
sec. 
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Quantitative Prediction of Earthquake Occurrence as 
Stochastic Phemomena 


By 
Keiiti Axt 
Geophysical Institute, Faculty of Science, Tokyo University, Tokyo, Japan. 


Abstract 


In this paper, the occurrence of earthquakes is regarded as a stochastic process and 
to it the theory of stochastic process is applied with an intention to investigate the possi- 
bility of the earthquake predictions by this procedure. At first the theory of prediction 
is roughly sketched. Discussions follow as to what quantity in the phenomena should be 
taken as the stochastic variable. As the result [n(t)—E(n(t))] is found to be the most 
appropriate quantity, where z(t) is the number of earthquakes belonging to specified range 
of magnitude which occurred ina specified area from t=0 to t=t and E(n(t)) is the average 
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of w(t) at t=t. 


The theory is applied to conspicuous earthquakes occurred near Shiriya- 


saki between 1900 and 1945, and the occurrence is predicted or extrapolated up to 1947 


and compared with the actual. 


Introduction—One of the essential problems 
in Seismology is to establish the law of oc- 
currence of earthquakes which, with the 
knowledge based on the present and past 
natural conditions concerned, can give us an 
information about their future occurrences. 
This law, different from the deterministic 
ones in various other sciences, must have a 
statistical character for the present, because 

‘of our present ignorance about the conditions 
by which the occurrence of earthquakes is 
influenced and the deterministic laws by 
which it is governed under such conditions. 
It is hoped that the statistical method used 
in the present paper will give us a certain 
macroscopic information about the phenomena 
which can not otherwise be obtained for the 
present. 

Many statistical studies have been made 
about the occurrence of earthquakes since 
the beginning of this century. The outline 
of these studies can be seen in the text-book 
by Matuzawa,“ for instance. In these studies, 
two sorts of method are important ; the one 
is to regard the occurrence of earthquakes 
as a Markov process as is done in the papers 
by WatanaBeE,® Kiracawa® and others. The 


other is to establish the periodicity in the 
sequence of occurrence of earthquake as is 
done by Matuzawa®). 

The method used in the present paper is 
somewhat different from the above two and 
is applicable to a stationary stochastic process. 
Therefore the occurrence of earthquakes to 
be studied by this method must be the one 
which can be regarded as stationary. But this 
is not always possible. For example, a series of 
after-shocks occurring exactly according to the 
Omort’s law can ‘ot be regarded as stationary. 
Another example is the so called trend, al- 
though this can be a part of a phenomenon 
which, if observed long enough, may well be 
a stationary one. In the former case, how- 
ever, we can regard the total series as an 
element of a stationary process if the time 
interval is taken long enough and if we give 
up to get any information about the process 
occurring within that time interval. In the 
latter case, we can pick up the stationary 
process only, by eliminating the trend from 
the total process. 

The material on which the studies will be 
made in the present paper is the series of 
210 conspicuous earthquakes which occurred 
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in the area of about 3x104km? near Shiriya- 
saki between 1900 and 1945. The data were 
supplied by the Central Meteorogical Obses- 
vatory, Tokyo®. 

Theory of prediction—Regarding a statio- 
nary stochastic process we have the follow- 
ing theorems), 

Theorem (1) A stationary stochastic process 
X(t), of which the mean value is zero, can 
be seperated into two parts as’ 


(1) 


where U(¢) has a continuous spectrum only 
and can be expressed as 


X(t) =U0)+ Vid), 


Uit)= Sb YVt—k). 


k=0 


(2) 


Here Y(t) is a stochastic process of which 
the mean value is zero and its every auto- 
correlation coefficient is zero. 

Vit) is a singular process which has a 
discontinuous spectrum only. 

Theorem (2) When there is no singular 
process in X(t), we have 


(3) X()= > by Y(¢—R) 
k=0 


(4) YQ)= za eXESR) 


where 


(5) Say Beug&0 (Rak 2: 35 we tee 


Q=oo=s1 and a= 


be=0 (R<O). 


The derivation of these theorems is in the 
papers about time series, e. g., OGAWARA’S®), 
Hereafter we shall consider the case when 
there is no singular process in X(¢), for the 
correlogram of our time series which will 
be given later has no definite singularity. 

It is important that the value of X‘¢) for 
the time ¢ can be expressed as a weighted 
sum of a series of Y(¢) for the time previous 
to ¢, where the mean values of X(¢) and Y(2) 
are both zero and every autocorrelation 
coefficient of Y(t) is zero, that is, 


E(Y(t)Y (t—k)) 
==) (Ral 2a 
DY) ( = 
Here’D*(Y )is the variance of Y(t) and E( ) 
is the operation to average the term in the 
bracket. 
Inserting (4) into (6), we have 


(7) DY)DX)> Qi Ox-i 
i=0 


(6) 


E(Y (#)X(t—k))= 


=={)) 
(R=172 30-6) 
E(X(t)X(t—k)) 
D*(X) : 
correlation coefficient of X(¢), and can be ex- 
pressed by the relation (3) as follows, 


where p.= This is the auto- © 


> Bidi+i 
(8) 


pr Sa 
Pas 

The relations (7) and (8) are important for 
obtaining the prediction operators a, and bh, 
from the known quantity p;. 

The predicted or expected value ae the 
time ¢+k (k=1, 2, 3, ....) is the conditional 
mean value of 

X(t+k) 
under the conditions 


X(t{—1)=a1_-4 Ose me se. 5 


where 2_; is the actual or sample value 
of the process for the time ¢—2. 


If 2; is known, y;-; can be obtained by the 
formula (4) as follows, 


(9) Vet ae tiick 
k=0 


Then the conditional average of x(¢+h) is 
given by the following formula, 

(10) EAX(t+k)) 
=2('S: bYE+h—1)+ 5 beste :) 


i=0 i=0 


Otherwise, from the formula (4), we have 
(11) lagi 


t=0 


=~ Sas E(XU+h-d)— ier e é 
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In this case the expectation is carried out 
step by step. (10) and (11) were called the 
extrapolation formulas of the first and second 
kinds in the work by Ocawara®, ~ 


The accuracy of prediction is measured by 
its variance 


(12) E[(X.¢+k)—EAXE+R)})] 
é 3 be DX Y) 
ui Bde 2D) 
where ‘ 
Om a 


In the practical application of the above 
theory, we obtain a, and dy not from the 
formulus (8) and (9) but by regarding the 
process of autoregressive scheme as follows. 
A stochastic process X(¢) of autoregressive 
scheme is one which satishes the following 
relation 


h 
(13) Sai X(t—i)= Vit), 
i=0 


where hi is finite, and Y(Z) is as before the 
process of which every autocorrelation coef- 
ficient is zero. Then we have 


h 

(14) PEL pra sles Bs er) 
Thus, when we have the correlogram of the 
time series, we choose a proper value h 
such that the assumed correlation function 
of this scheme of rank h chosen may 
interpret the whole correlogram, and solve 
-the system of initial # equations in (14) and 
obtain az. And then bd; is calculated from 
(5). 

Refering to the above theory, Woop“ dealed 
with the time series in a similar way as 
above, and also WIENER) referred to the dis- 
crete time series in his book ‘‘Cybernetics”. 

Application.—Before the application of the 
above theory to our problem, it must be 
considered what quantities in the phenomena 
should be taken as the stochastic variable. 
The idea of taking the energy of earthquakes 
as the variable seems reasonable because it 


has an obvious physical meaning. But the 
energies of large and small earthquakes dif- 
fer considerably and the largest earthquakes 
play a predominant role. It is therefore 
dificult to regard the series of energy as 
stationary during our period of observation. 
For the present, therefore, more preferable 
one is the number of earthquakes occurring 
in a specified area in a specified range of 
magnitude. In the present paper, the follow- 
ing three expressions based on this number 
are discussed and the above theory is applied 
to the best of them: 


(a) series of time interval between suc- 
cessive earthquakes. 
(b) series of number of earthquakes per 


unit time e.g. per month. 

(c) n(t)—E (n(Z)) 
where n(¢) is the number of earthquakes 
between ¢=0 and ¢=/ and E(n(Z)) is the 
average of mnt) at t=¢. If a@ is the 
average number of earthquakes per 
unit time, H(n(t)) is equal to at. 

In the followings, these three time series 
will be constructed from the data referred to 
before and their statistical property and the 
applicability of the theory of prediction will 
be examined. 

(a) Series of time interval between suc- 
cessive earthquakes. 

The principal probability density of this time 
series is shown in Fig. 1. It seems to be in 
a good accordance with the exponential 
density curve also shown in the figure. But, 
for the effective use of our method, it is re- 
quired that the distribution function should 
be normal i. e. Gaussian. Therefore our 
series must be transformed into a new one 
of which the distribution function is normal, 
and the correlogram of this new series must 
be formed. The result is given in Fig. 2. 
The figure shows that in the correlogram 
for the first period (1900-1910), the hypothesis 
that all the autocorrelation coefficients are 
zero cannot be rejected on the 10% level 
of significance except one which may be 
assumed as positive. This is the same in the 
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second period (1910-1920), but in this case 
one significant autocorrelation coefficient is | 

0.5 negative. It is almost hopeless that we can, | 
from such a sample correlogram, get any 
theoretical model of autocorrelation function 
for a non purely random process to which 
we shall apply the theory. So this expression 
of occurrence of earthquakes should be ab- 
andoned. 

(b) Series of number of earthquakes per 
month. The principal distribution of this 
series is shown by dots in Fig. 3, in which 
the curve shows the Poisson’s distribution of 


0 5 10 —> Months mean value 0.4. In this c se the correlogram ° 
Fig. 1. Principal probability density of is formed without the normalization of dis- 
the series of time interval tribution. Since the number of earthquakes 


per month is almost always one or zero, the 
correlation coefficient of those numbers 
which are mutually ¢ months apart has an 
[ 5 To almost linear relation to the probability that 

oe an earthquake occurs ¢ months after the time 
ay 1901~1910 when an earthquake took place, i. e. the sim- 
1 plest transition probability, and we can test 
if me the uniformity of the probability about months 


aie Sony eek Srp al by the method of % square test and at the 
c= ia me 

TF same time the significance of the correlation 

ee coefficient. The result is that the correlation 

= coefficients are more significant than in the 


Fig. 2. Correlograms of the series of time focmer case. But as shown in Fig. 4, if we 
interval 1 


—> Months 


1901~ 1910 


0.5 0 


30 


1911~1920 


1921~1930 
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Fig. 3. Principal distribution of the series Fig. 4. Correlograms of the series of number 
of number of earthquakes per month of earthquakes per month 
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Fig. 5. 2(t)—H(n(t)) obtained from the data which contains 210 major earthquakes 
oecurred near Shiriya-saki from 1901 to 1945. 


fit continuous autocorrelation functions to 
those sample correlograms, the sample points 
must vary largely about the function. And 
, other more suitable expressions are required. 

(c) n(t)—E(n(t)). The meaning of this 
quantity was given before. The actual vari- 
ation of this function for our sample is shown 
in Fig. 5. For the convenience of analysis 
this is seperated into two parts; one is along 
period variation and the other is a short 
period variation. Smoothing the original 
variation by the method of moving average, 
the former is obtained and subtracting this 
from the original the latter is given. These 
two variations can be regarded to have a 
negligible correlation, and the seperation is 
justified. 

(1) Short period variation. The principal 
probability density of this variation is given 
in Fig. 6 and seems in a good accordance 
with the normal probability density of mean 
value 0 and standard deviation 1.4. The 

“correlogram can be immediately formed in 


Fig. 6. 


Principal probability density of the 
short period variation 


this case and is show in Fig. 7. In this case 
the correlograms.are sufficiently smooth that 
we can reduce significant and continuous 
autocorrelation functions from them and 
regard them as the samples of these functions. 
But, as shown in the figure, these four cor- 
relograms are considerably different from 
each other and it is not preferable to regard 


—e Months 


1901~1910 


1911~1920 


1936~1945 


Fig. 7. Correlograms of the short period 


variation. 

them as samples of a single autocorrelation 
function. The correlogram for 1910-1920 
has no negative values and diminishes 
almost exponentially, and others are like 
a damped oscillation. From a very simple 
idea that the correlogram for 1921-1936 and 
that for 1936-1945 are comparatively similar 
and that the correlogram after 1945 will 
be the similar one, the prediction operators 
are obtained from the correlogram for 1921- 
1936. 

Assuming that the correlogram is one of 
the autoregressive scheme of h=3, and tak- 
ing unit time of five months, we get from 
the correlogram 


then we obtain from the formula (14) 
a=—0.18 a,=0.42 a3=0.12 
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and from the formura (5) 
b,=0.18 b,==—0.39 b3=—0.27 


The standard deviation which gives the error 
of prediction is for the first step prediction 


(i= NAVE 
and for the second and third steps 
fR=T.10 Jp heed, 


The result of prediction for the shert period 
ariation is shown in Fig. 8. 


predicted -curve 


Fig. 8. Precicted curve of the short 
period variation. 


(2) Long period variation. The correlo- 
gram of the long period variation is given in 
Fig. 9. If wetake the unit time of five mon- 
ths as inthe case of short time variation, we 


Fig. 9. Correlogram of the long period 


variation. 


can hardly get an autocorrelation function 
of autoregressive scheme of relatively small 
h which should interprete the whole correlo- 
gram. In order to avoid this difficulty a 
method developed by Ogawara is used. At 
first letting the unit time be a necessary 


length +, we solve the next equations and 
obtain a’; , 
h 
COpwee eee ni=0 (h=1,2)3, 20,8) 
v1=0 


from the definition of the autoregressive 
scheme. The next equations hold for k>h-+1 


h 
(16) Das pr-+=0 : 


The characteristic equation corresponding to 
this difference equation is 
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(17) 


ls ie 
Sage" *=0 
i= 


Let h roots of this equation be 


; #2nif 4 
BossPie= 272... -. > By, 1e 4+—* for odd h 
2 


By et2nis1 ; By et2tit 2 —— 


for even h, 


where f/f; is the component frequency and 
is the component damping ratio contained 
in the correlogram. Then we set up a h-th 


order algebraic equation of which hk roots. 


sir 
are, e. g. for odd h, Bs, Pr geet 


My se eng 
n=l 

er) grr, 
n 2 


a (Ale Betyehent) 


Br-1 2tf;, 1 Bn-1 =2nif;, 1 
n 


Se: ae 27-e@ 3 
n 
=2%+4,2°-14 |... +an_12+4p 
And this gives us a which is the prediction 
operator of unit time of ct/z. If it happens 
that By is negative in the case of odd h, we 
introduce two roots —fye*** in place of fy 
and the rank h becomes f+1. 
Following the above method, we take rt of 
0.4 years and get from the correlogram 


01=0.1 0'2=—0.3 o’3= —0.3 
then we obtain from the formula (14) 
a,=—0.01 a’,=0.28 a’3=0.25 


and the prediction operator of the unit time 
of ten months is given as follows, 

@=—3.41, a,.=4.68, a3=—2.92, as=0275 
The error of prediction, i.e. the standard 
deviation of expected value, for the first step 
is 2.1, and that for the second is 3.1. 

The result of prediction for the long period 
variation is shown in fig. 10. The curve in 
the figure shows an unexpected steep des- 
cent, for n(¢) is from its definition a mono- 
tonous increasing function and the gradient 
of n(t)—E(n(¢)) must be larger than —@w where 
a is such that E(n(¢))=at. But here ignor- 
ing such somewhat unreasonable fact that 
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1915 


Fig. 10. Predicted curve of the long period 
variation hs 


the gradient of predicted curve is less than 
—a, the predicted cuves of long and short 


‘period variation is put together and shown 


in Fig. 11. And for comparison, the actual 
curve is given in the same figure. 


predicted curve 


actual curve 


Jan. 191A 


Bigs Ti. 
Remarks—The whole amount of informa- 
tion obtainable by this method about the 
phenomena, which is got as stated in the 
introduction with all our ignorance about 


Predicted curve of x(t)— E(x(t)) 


* conditions and laws unknown to us for the 


present, was given in the correlograms so 
far obtained. And the prediction was based 
on these correlograms. Though the predi- 


-¢tion failed partly in the present case, the 


writer hopes that such a quantitative method 
for the study of occurrence of earthquake 
becomes a clue for a further approach to 
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the problem. As far as the present case is 
concerned, the accuracy of prediction is 
damaged mainly by the relatively large 
variance of the long period variation, and 
the ambiguity in the choice of a correlogram 
from four ones in the case of short period 
variation. These difficulties may be over- 
come step by step in the future by con- 
sidering the correlations which may exist 
between the occurrence of earthquake and 
other related observable phenomena. 
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